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I. Introduction 


In remote sensing of the Earth with microwaves^ we are inter- 
ested in surface and subsurface geophysical features rather than 
atmospheric effects. Hesearches are being performed, on ocean sur- 
faces^"^^ as well as in land areas. The articles by Staelin^^ 
and Tomiyasu^*^ provide a comprehensive review on atmospheric and 
surface observations, passive and active sensing, and their impli- 
cations and applications. 

Our interest is directed toward land areas? especially snow and 
ice covered land or water, and deserts. The objective is to develop 
theoretical models that are practical and useful in data interpreta- 
tion . 

We now briefly review the relevant theories that are of direct 

48 

interest to us, Gurvich, Kalinin, and Matveyev derived expressions 
for the brightness temperature of a half-space random medium v;ith 
laminar structure and used to interpret Kosmos 243 data for Antarctic 
areas. Assuming a half-space containing random distributed scatterers 
embedded in a low loss dielectric, England* examined emission dark- 
ening caused by scatterers with a radiative transfer approach. 

Stogryn examined the brightness temperature of a vertically struc- 
tured medium with no scattering. He also studied scattering by 
random dielectric constant fluctuations in the low frequency limit 
using a perturbation approach. 
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II. Accomplishments of the Past and Proposed Work for the Future 

The dominant, factors that determined remote sensing measure- 
ments are absorption, layering, anisotropy, submerged scattering 
centers, nonuniform temperature distributions, and suirface roughness. 
Assuming uniform temperature distribution, the first three factors 
are accounted for in tlie composite model, v;hich also extends the 

half-space model of Gurvich et al to treat layered media. Tlie com- 

53-54 

posite model was used to match (Appendices 1, 2) spacecraft 
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data with ground truth measurements. Satellite data were also 
studied with a stratified model. (Appendix 3) Gurvich et al as- 
sumed uniform temperature distributions. We solve the problem of 
microwave thermal emission from a half-space random medium with 
nonuniform temperature distributions. (Appendix 4) The observed 
data for slab medium exhibit coherent effects. We developed modi- 
fied radiative transfer equations and solved for the reflectivities 

57 

of a -two layer random medium. (Appendix 5) The results are ap- 
plicable to active sensing as well as passive sensing. For a 
stratified medium with a non uniform distribution, we use the fluc- 
tuation-dissipation theorem and solved for the brightness tempera- 
5 8 

ture. (Appendix 6) All these results are easily used. They are 
summarized in tiie subsequent sections. 

The task is, however, by no means completed. We have concen- 
trated on structure.*? that are laminar. Three dimensional scattering 
needs to be considered. Surface roughness has so far been neglected. 
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Its effect should be included. Our final goal is a composite model 
which accounts for all the determining factors. We would like to 
be able to identify the effect of each individual factor on the 
observed data. The theoretical models should be used to interpret 
all available data. Only through data match can a useful and prac- 
tical theory be developed. 




III, stratified Media V7ith Uniform Temperature Pi s tr ib utio n 


The emissivity e is related to the reflectivity r by 

♦ 

e = 1 - r, (1) 

The relfectivity is equal to the square of the reflection coefficient 


r = 1 - . 



( 2 ) 


For the stratified medium as shown in Figure 1, R is calculated 
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in the form of continuous fractions 
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whete for horizontal polarization 
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The above result is exact, Absorption is taken into account by 
the imaginary parts of y and e. 
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Figure 1. Stratified medium with uniform temperature distribution 






Half-Space Random Media with Nonuniform Tempearature Distribution 


For a half-space medium with constant absorption and scattering 

coefficients, the brightness temperature is determined by the fol- 
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lowing simple closed- form formula 


(a t - r^^lot - K^) 




where we assume that the temperature distribution has the profile 




In (6) r^^ is the reflectivity of the medium, 0 the observation 

angle, the absorption coefficient, and 
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In (8) ic is the scattering coefficient, p, the forward- 

o X 

scattering fvinction, Pj^ the backward-scattering function, and the 
scattering albedo is 




For a medium with no scattering, ic^ - u « 0 * Thus 
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For a medium with scattering and uniform temperature distribution, 
T, ~ 0, Thus, 
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( 11 ) 


For uniform temperature distribution and small scattering albedo, 
we expand (11) to first order in w and obtain 


T^ = {1 - r„. ) 1 “ 
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cos 0(1 + 4kj^^Jl^cos^e) 


Niimerical results are given in Fig. 2 to compare brightness tempera- 
ture for TE and TM waves as a function of radiometer viewing 
angle 0^, witli tenperature profile T(z) = 222 H- 34 exp[0.81z]. 
Clearly, scattering lowers the brightness temperature at all viewing 
angles. • . 
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Figure 2 
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Two Layer Random Media with Uniform Temperature Distribution 


Consider a slab of random medium of thickness d. A plane wave 
is normally incident upon the slab. The reflectivity r is solved 
from tlie modified radiative transfer (MRT) equations, which include 
coherent effects due to the two boundaries 
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For a half-*space random madium, we let d ^ co. 


= ^ 01 + 
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which is identical to previous results. When the second layer is 
a perfect conductor, - -1, r ^2 “ ^ have 
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It is interesting to note that as u + 1, we have a *>• 0, f ^ 1 . 

Thus r 1 and all the incident power is reflected. In Fig. 3 we 
illustrate the emissivity for a slab random medium. We see that 


IZ 


Figure 3 . 

Smissivity s= 3,0(1 + i 
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■a 

a = 2 mm,. 
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scat-tering dampens the interference pattern and decreases emissivity 

in general. The existence of the interference patterns depends on 

the location of the subsurface and the extinction loss of tlie random 

itvadium, ... 

The brightness temperature of a stratified medium with nonuni- 

5 8 

form temperature distribution is calculated in the following form 
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VI , stratified Media with Monmiform Temperature Distribution 


The brightness temperature of a stratified medium with nonuniform 
temperature distribution is calculated in the following form^® 
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Here is the observation angle, e^" is the imaginary part of the 
permittivity in layer t. Tp(z‘) is the temperature distribution * 
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in layer 

region 

matrices 


and and are wave amplitudes in 

that are related to those in other regions by propagation 
59 
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A tv?o layer composite model'* is used to match micro- 
wave emission data of fresh water ice and snow in order both 
to test the model and to interpret the data. The magnitude 
of coherent effects are empirically estimated with the object 
of explaining features of the emission spectra, and the model 
is modified to deal with polarization effects and oblique 
incidence . 
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I . INTRODUCTION 


Remote sensing of the earth at microwave frequencies is 
of great use because of its all v/eather and day-night opera- 
tional capability. This is due to high transmission through 
the atmosphere and cloud. Many factors such as absorption, 
layering, scattering, polarization, incident angle, anisotropy 
and medium inhomogeneity affect the radiobrightness spectra. 

In order to interpret the data, all these effects should be 
incorporated in a theoretical composite model. Development 
of such models, combined with the interpretation of the measure- 
ments V7ill lead to numerous practical applications. For 
remote sensing of ice and snow, applications include ice 
dynamics, climatology, glaciology and navigation etc. 

Effects of scattering due to medium inhomogeneity have 
been accounted for by A.S.Gurvich et all^^ The theory enables 
us to distinguish the spectra of the continental glaciers frpm 
those of shelf glaciers. Inversion of the data is performed. 

The estimated values for the variance of the density fluctu- 
ations and the characteristic length of ice crusts agree with 
the measurements. In order to develop a more complete model, 
the formulas given in that paper have been extended and com- 
bined by J.A.Kong with the stratified model which also 
accounts for the effects of layering, absorption and aniso- 
tropy. 
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In the present paper, the model developed by Kong is 

[ 3 ] 

referred to as the composite model and xs used extensively 
to match emission data* The problems encountered in the 
process of matching and inverting the data suggest the pres- 
ence of additional factors that influence emission properties. 
With the object of improving the model, the magnitudes of 
coherent effects are estimated and the tv/o layer composite 
model is modified to deal with polarization effects and 
oblique incidence. 
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II . DATA MATCHING 


The Composite Model 

f31 

The two layer composite model consists of a layer 
of snow with depth d and relative complex permitivity ~ 
e£ +ie£ on top of ice, water or land. 

air 

d snow 2 q, = ei +iej^ 

T ice, water or land ^2 £2 ■*'^^2 

Pig. 1 Two layer model 

The reflectivity and scattering factors are calculated sepa- 
rately to give total emissivity e of the structure: 

e = 1 - r - s (1) 

where r is the reflectivity of the object and s is the power 
coefficient of scattering effects. Observations are assumed 
to be made from nadir. The oscillatory behavior of r and 
scattering effects due to the substrata are neglected in the 
model. The emissivity depends on frequency and seven model 
parameters; d, € 2 ' ^o ^1 ^2 

the relative complex permitivity of a top and a bottom layer 
respectively, is a characteristic correlation depth for 
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the refractive index, P is a dimensionless scattering para- 
meter and is related to the variance of refractive index 

1^1 * o 

fluctuation by P - — C/no' where n| ***nji = 

* 1/2 ^l" ^ 

Data Matching 

The model has been tested on the emission data from 

F41 

fresh water ice and snow published by T, Schmugge et al . 
The data was taken by microwave radiometers in the frequency 
range 1.43Ghz to 37.5Ghz during the flights made over dry 
snow with different substrata: lake ice, wet soil, glacier 

ice, firn and wet snow. The lowest four of six frequency 
channels are the nadir views. The 37.5Ghz channel provides 
vertical and horizontal polarizations with the look angle 
of 45®, The information on the second highest frequency 
channel {19.35Ghz) is ambiguous, but the nadir observation 
is assumed for our purpose of inverting the data, Emissivity 
as a function of frequency is calculated for each set of 
different parameter values in order to fit the theoretical 
curve to the data within the margin of 0.03 in emissivity. 

We have assumed the constant depth and ground temperature. 

The calculation is made using the computer language, I^ACSYMA? 

^Project mac’s symbolic MAnipulation system; developed by 
the Mathlab Group, Project MAC, MIT, whose work is spon- 
sored by the Advanced Research Project Agency, Department of 
Defense, under contract number N00014-70-A-0362-0001 
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The data and the results of inversion are presented in 
Appendix I, and explained for each site individually in the 
following paragraphs. • 

i) Bear Lake ( d == 0,15m, altitude 1805m) 

[model] [measurement] 

0 I air 


r- 

1 0,15 snow 0.15 



ice 0.25 


i 

water 

Fig. 2 Model and ground Measurements 

The data shoves a sharp drop in emissivity at low fre- 
quencies, This characteristic is best matched by setting 
c£ == 3,2 and €2 = 80 + 80i, For such dielectric properties, 
the penetration depth at l,43Ghz is about 47,5cm, The sharp 
drop can thus be accounted for by layering, especially due 
to the large values of ^2* This is in good agree- 

ment with the ground measurement which indicates that there 
is water beneath snow. The modeling is crude because v?e 
have ignored a slab of ice between snow and water. The para 
meter values obtained for snow are very much like those of 
ice. We believe that the model should be extended to a 
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three layer model in order to account for the dependence of 
the data on both snow and ice layers and scattering effects 
caused by ice, . 

ii) Steamboat Spring ( d ^ 0.8m, altitude 2070m ) 

[model] [measurement] 

0 ^ air 

1 0.8 snow 0,8 

2 f land 

Fig. 3 Model and ground measurement 

As mentioned in the Bear Lake case, low emissivity in 
the data at the lovrest frequency channel (1.43Ghz) is caused 

■k 

by the substratum. Penetration and probing depth at 1.43 

gigahertz are calculated to be 39.4cra and 5.3m respectively, 

while d is o,8m. The values 60-80 and 0-80 for ej ^2 

suggest that v/et soil with moisture of greater than 35% 

appears to be water-like from the radiometers. The value 

= 2.2 + 0.04i is typical for snow. This agrees with the 

measurement. It should be noted that Z is 0.4mm, the lowest 

o 

among all six tested. It indicates this snow has fine struc- 
ture. 

'"‘''The maximum depth that we can probe depends on sensitivity 
of radiometer. Assuming the sensitivity of 2°K, we define 
probing depth as the maximum that can distinguish a two layer 
medium from a half space case. 
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1^ 

iii) South Cascade La)ce { d = 0*5m, altitude 1610m ) 



tmodel] 

[measurement] 



air 


0 

_ L 

light dry snow 


1 

0,5 

light snow with 

! 

liquid 0.3 

2 

t 

wet snow, slush. 

I ^ 

water 

Fig. i 

Model and ground measurement 

• 


At this site, the simulation of the real structure is 
very crude. The largest values for the imaginary part of 
complex permitivity of snow are obtained. This reflects 
the presence of liquid vrithin snow as indicated by the ground 
measurement, V7e also note that there is very little volume 
scattering (P - 0.35) compared with other sites. This is 
in agreement with Schmugge’s explanation that scattering 
effects are washed out by the increased effective loss tangent 
for snow when liquid water is present. 


iv) South Cascade Glaciers 


P-0, 

d = 

4.9m, 

altitude 

1610m, 

below 

firn 

line 

P-1, 

d - 

6.8m, 

altitude 

1890m, 

above 

firn 

line 

P-3, 

' d * 

8.4m, 

altitude 

2040m. 
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[model] [measurement] 



Pig. 5 Model and ground measurement 

The parameter values obtained indicate that there is 
no liquid content in snow. The effects of volume scattering 
for all three places at this site are relatively large (P = 

f5i 

1.8-2. 3), and this checks with the past observation •• that 
multi-year ice has lower emissivity at higher frequency chan- 
nels than new ice. The values for ^2 distinc- 

tively different from those of other three sites and this is 
justified by the fact that there is ice or firn, instead of 
water or moist, land underneath the snow. Because we have 
ignored the scattering effect caused by the substratum, the 
model can be modified to deal with such effects. 
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Summary 


The parameter values shown in Appendix X are obtained , 
as a result of tedious and time-consuming computer manipu- 
lation to fit the data within the margin, Emissivities at 
different frequencies are calculated for a set of parameters 
and compared with the data by inspection. This procedure 
is repeated unsystematically until we find a reasonable 
match. Consistent results cannot be guaranteed because 
such a method often depends on the person doing the manipu- 
lation. D.Chu has wrtten a computer program^ ^ to perform 
such inversions automatically. The program has been tested 
on the same data. Although it is not as flexible as a 
human operator (the program approaches the data only from 
above and belov?) , his method has advantage of providing con- 
sistent results in a short time. 

The composite model proved to be capable of explaining 
many features of the emission data. These characteristics 
can be related to physical parameters by matching the data. 

The results of the inversion suggest that we are able to deter- 
mine the dielectric properties of a top and a bottom layer, 
relative magnitude of liquid content, size of the scattering 
centers and relative magnitude of scattering effect. We 
have not tried to invert the thickness, but we believe that 
we can estimate d in case of ice or snow, provided that d is 
smaller than 0.5m. 
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The above mentioned parameter values, for example, can 
distinguish dry snow from wet snow, first-year ice from multi- 
year ice, snow or ice with fine structure from that with 
rough structure and of course, snow and ice from water or 
land. Although we dealt only with smooth surfaces, the effect 

of surface roughness can be explained by modifying the reflec- 

[71 

tion coefficient used in the model 

The major difficulty in matching Schmugge's data was 
caused by his ambiguity about the incident angle of the 19.35 
gigahertz channel. This channel provides the horizontal 
polarization and scanning mode. The abrupt drop in emissivity 
at this frequency made the matching very difficult. if we 
assume that the data was obtained by scanning the radiometer 
from 0® to 45® for example, this sharp drop can be accounted 
for by the polarization and angle dependence of emissivity. 

The present model is limited to only two layer media. 

For some data, more sophisticated modeling is necessary. 

■In order to explain the data from polarization channels, the 
model should be modified. This is done in the following 
part of the paper together v/ith the examination of oscillatory 
behavior of the reflectivity. 
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II, MODIFICATION OF THE MODEL 


Folarigation and Oblique incidence 

Dual polarization and scanning capability of microwave 
sensors are attractive features which can provide additional 
information on radio emission. Two of six radiometers used 
for Schinugge ' s data have polarization channels : dual at 

37,5Ghz and horizontal polarization at 19,35Ghz. The present 
model, however, is insufficient to account for the measure- 
ments on these channels. In order to deal with effects on 
emissivity of scanning and polarization, we modify the com- 
posite model as follows. 

We assume the angle of incidence and the angle of 
refraction as shown in Pig. 6, 


X 



Pig. 6 Two Layer Model (Oblique Incidence) 




olbiLOT 0 ? 
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Prom the composi€e model, reflectivity r is given by 


Poif 


-4klitd 


r = 


( 2 ) 


1 ^|K0lP-Pl2f- 


0 1/2 

where + k^^” i = (e| + ej i - Sin 9^) , 

and Rp^ and ^^2 reflection coefficients of a top and 

a bottom interfaces. For the horizontal polarization. 


1 - 


klx 


koK 


ROI = 


1 •!■ 


kix 


and 


1 


kox 


l-2x 


klx 


1^12 - 


1 t 


k2x 


^Ix 


v;here k 2 jj = (e^ t £2 i " Sin 0^) 

^ox ^ “ Sin^e^)^'^^. 

For the vertical polarization. 


1/2 


(3. a) 


(3,b) 


llOl = 


eik 


ox 


- k 


lx 


5lkox *** kix 


(4. a) 


R12 = 


^2kix “ £lk2x 
^2kxx + ^lk2x 


(4.b) 


S£.'SS“- 
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The polarization dependent power coefficient of the 

scattering effect due to medium inhomogeneity for an unbounded 

n 1 

half space is given by Gurvich et al . By extending their 
result to a two layer case, we obtain. 




n 


2n? 


(T, 


no 


1 + (2n|kjjZ^Cosej^) 


Cos 0, 


X (1 -|Rq/) • ( 1 


-4 'koni- d/Cos©i 
0 / 


(5) 


where Sj^ is the scattering factor for the horizontal polari 
z at ion, and 


®v ^ ®h 


( 6 ) 


for the vertical polarization v;here 0^ and are related by 
Snell's law: 

nj[*Sin - Sin (7) 

The emissivity for each polarization are, 

e, « 1 - r. “ s, (8. a) 

h h h 

and 


e = 1 - r “ s (8,b) 

V V V 

where reflectivities r and r„ are calculated for polarize- 

h V 

tion from (3) and (4) respectively. 


44 


We have used the parameter values matched before 
(Appendix I) to calculate the emissivities at 37,5Ghz chan- 
nel for both polarizations. Appendix II lists the results 
in comparison with the actual data. in Pig. 7, we show 
the incident angle dependence of emissivity for the two 
polarizations. 



1% 


§ 

't 

Coherent Effects 

The composite model assumes all radiation is rendered 
incohernt by scattering and medium inhomogeneity,- In general 
scattering and absorption attenuate the coherent wave. In 

■j: 

the case of ice or snow, we observe a large scattering effect 
in the higher frequency channels. Absorption is relatively 

■ l 

small for a thin layer medium and at low frequencies. f 

In such situations, we expect the data to show some charac- 

teristics of coherent radiation in low frequency channels. ■ 

Since such characteristics cannot be matched by the present ’ 

model, we are to remove the assumption of total incoherence 

from the model in order to come up with better match and to 

investigate the magnitude of coherent effects. 

The exact expression of the reflection coefficient for :■ 

a two layer medium is, 

H. ” — - - — 

where R,., , R . k. are as defined earlier (Fig. 1) . 

01 12 lx 

By separating the real and imaginary part of ws 

0 

the reflectivity r =}R| : 

2N Cos2kiid 

1 2N •lRo]j-iRi2l2"^^l^'^- Cos2ki^^d 


(9) 


find 


(10) [j; 
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where N characterizes the magnitude of the oscillatory 
behavior of the reflectivity. If we assume the wave to be 
incoherent/. N is 0, A totally coherent wave in this case 
corresponds to N - 1. By adjusting the coherence factor 
N, we have matched the data of South Cascade Lake (Appen- 
dix I). We have obtained the best result with N = 0.25. 

The imaginary part of accounts for absorption of coher- 
ent wave. The scattering loss is neglected because coher- 
ence effects are important only at the lox'/est two frequency 
channels where there is negligeble scattering. The match 
is improved and the graph is shown in Pig. 8 in comparison 
with our result before. Pig, 9 shows the emissivities with 
different values of N. 

Although the above procedure accounts for some features 
in the emission data/ the method is purely empirical. The 
theory that accomodates coherent effects has been investi- 

r 81 

gated by L. Tsang . In his paper incoherent wave and 
scattering effects are combined with a coherent wave which 
suffers scattering loss as well as absorption. 
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Emissiviijy 


€- = 3*2 + 0e06l 

ej = 80 
0.8m 
P ^1.28 
= OoSmin 




/ / / 

/ / 

/ / 


liao B«aoB # cj o o « 


curve without oscillation 
envelop of oscillation with 
envelon of oscillation with 


Frequency 


Fig, 9 Emissivities with T>ifferent Coherence Factors 
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IV. SUM14AHY AND SOGGESTIONS FOR THE FUTURE STUDIES 


We have found that the two layer composite model can 
provide useful parameter values from emission data. The 
inversion is a tedious process. The theoretical curve 
must be matched with the data as closely as possible in 
order to obtain the reasonable values. The modification 
to allow some coherent effects in the model has improved the 
matching. The generalised model which deals with polari- 
zation and angle dependence accounts for the additional 
features of the emission spectra. We need more accurate 
data, however, in order to perform data inversion using the 
modified model. Ambiguity about the oscillatory behavior 
can only be cleared up by probing at many different fre- 
quencies. Also if the scanning angle had been given on the 
19.35Ghz channel of Schmugge's measurements, the data match- 
ing would have been far more succesful. Since polarization 
can provide additional information, it is desirable to equip 
v;ith dual polarization for all channels . Such a feature is 
particularly useful when we study surface areas. 

It is essential to compare the values obtained using 
the model with the actual ground measurements. Except for 
the snow or ice thickness, however, Schmugge's measurements 
of geophysical parameters are ambiguous. In order to 
further the study .that relates the model's parameters to 
geophysical parameters, precise ground measurements are 
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required. Small scale experiments on the detection of radio 
emission can be set up for such a purpose. 

The model can now deal with absorption, anisotropy, 
layering, polarization, incident angle, scattering and medium 
inhomogeneity. As for layering, we have only tried a tv/o 
layer case. Often, a two layer model becomes a crude approx- 
imation of a real situation. Because the present model 
neglects scattering due to the underlying medium, a three 
layer model should replace the two layer one when the substra- 
tum is ice and when the whole structure is more complicated. 

It should be pointed out that data are expressed in emissivity 
which assumes the constant ground temperature. This approxi- 
mation is valid when the length of temperature wave is consid- 
erably small compared to the thickness of the layer. In 
this case temperature variations caused by vreather changes 
in a thin surface layer have little effect on emissivity. 

For the more general case, the effect of temperature variation 
vrith depth can be studied. 
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APPENDIX I 



[Results of 

the Inversion] 



Bear Lake 

Steamboat 

S. C. Lake 

d(m) 

0,15 

0.8 

0.5 


3.2 

2.2 

2. 9-3.2 

ej 

0.04-0.048 

0.04 

0.095-0.11 

^2 

80 

60-80 

70-80 

^2 

80 

0"*8 0 

0-20 

P 

1.1-1. 6 

1.34-1.38 

0.35-0.45 

z (iran) 
o 

0.8-1. 2 

0.37-0.4 

0.5-0.65 


S.C.G,P-3 

S.C,G,P-1 

S.C.G,P“0 

d(m) 

8.4 

6.8 

4.9 

4 

3.1 

2.7-2. 8 

2. 6-2.7 


0-0.0004 

0-0.0004 

0-0.0008 

^2 

3. 0-3.2 

3. 1-3. 2 

3. 0-3. 2 

e'* 

0.004-0.05 

0.004-0.05 

0.03-0.05 

P 

1. 8-2.1 

1. 9-2.0 

1.75-2.3 

Zq ( mm) 

0.75-0.9 

0.9-0.95 

0.5-0. 8 

Thickness 

d is given by the measurements. 


The data 

and the theoretical 

curve are plotted 

in the follow 


ing pages. 
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[Emissivities at 37.3GIIz for both polarizations] 

\ 


Actual data 

observed Figures calculated 


Location 

Horizontal 

Vertical 

Horizontal 

Vertical 

South Cascade 
Lake 

,856 

.941 

.830 

.964 

Bear Lake 

.739 

.862 

.767 

.928 

Cteamboat 

Spring 

,793 

.862 

.800 

.953 

Glacier 

P-0 

,768 

.867 

.751 

,924 

Glacier 

P-1 

.799 

.888 

.767 

.929 

Glacier 

P-3 

.769 

.877 

.731 

.906 


Radiometers are aimed at 45° from nadir. 

For the calculation, we have used the results of data 
matching given in Appendix I, 

eeprobucibility 
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ABSTRACT 

A computer program is described which interprets data 
according to a mathematical model. The program begins 
with a wide range for each of the unknown model parameters 
and tightens that range about the actual value of the 
parameter. Examples are presented which demonstrate the 
successes and the fail'ings of the program. To show the 
utility of the program, emissivity measurements, published 
by T. Schmugge, are inverted according to J.A, Kong's 
Composite Model. The results agree with those obtained 
after much effort by M. Nagase. 
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INVERSION 


r 


A mathematical model relates an observed quantity 
with a set of parameters. 

observation = f (parameter^, •••, parameter^) 

If the value of each parameter were known, the observed 
quantity could be found. In remote sensing, however, 
observations are made a different values of one parameter, 
such as frequency. Then the remaining parameters are de- 
termined from the observations. Operating a function 
in reverse this way is called inversion* and involves the 
solution of (n-1) Simula taneous equations. Because these 
implicit equations can be complicated, inversion may be 
difficult . 

The Composite Model was devised by J,A, Kong to in- 
terpret the thermal radiation from masses of snow and ice 

in terms of structural and dielectric properties.^ 

« 

Emlssivity, the measured quantity, is modeled as a function 
of five parameters and frequency, 

emlssivity = f( e ’ , e”, p, z^, d, 

Prom measurements at five frequencies, these parameters 
can, in principle, be determined. Exponential and poly- 
nomial expressions in the model, however, make the solution 
prohibitive, 

ande" are the real and imaginary parts of the dielectric 
constant, p aud z^relate to the magnitude and frequency 
distribution of scattering. d is the thickness of the ice 
or snow, and tois the angular frequency. 
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It would be convenient to have a computer program to 
perform the inversion. The computer would take any model 
and a set of data and return values for the model para- 
meters. Such a program would be particularly useful for 
evaluating models. If the parameters had physical signl- 
ficance, they could be measured independently and compared 
with the inverted results. A program is described here 
which uses an iterative process to take a wide range of 
parameter values and tighten it about an inverted value. 


PRINCIPLES and PROBLEMS 


The further a parameter is from its true value, the 
worse the prediction will fit the data. This is the prin- 
ciple behind the program. ■ Two sets of limits, and 
are chosen to generate curves, one above and one below 
the nieasurements . As the range of each parameter is de- 
creased, the limiting curves approach the measured curve. 

By minimizing the difference between the. prediction and the 

measurements, the program tries to bring each parameter 

♦ 

closer to its proper value. 



There are two problems with this approach. Pirst, 
a parameter may not yeild the minimum error at its true 


value unless all the ot* er parameters assume their actual 
values. Second, if changing either of two parameters has 
the same effect, how can one decide which parameter to 



If the model were gerfect, when all the parameters 
took on their true values, the variance between the pre- 
diction and data "ould be aero. Then, i^ one parameter, 
for instance e’, were allowed to change, the minimum in 
the error would increase. Were a second parameter allowed 
to change, the minimum in the error curve would, in general, 
be shifted from the actual parameter value. 



During the inversion, none of the five parameter values 

- ■ y 

is known. The error gives only a crude idea of where the 
parameter lies. Rather than declare the value of e’ to 
be four, the program uses the error to decide which of the 
limits is nearer the true value of e'. As the parameters 
converge to their proper values, the minimum in the error 
curve should approach the true value of e*. 

The effects of the parameters generally overlap. For 
instance, decreasing the emissivity at l0.7 gigahertz 
'may be'done by increasing the magnitude of the scattering 
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REPRODUCIBILITY- OF THE 

original page is poor 


The first choice 


2 

curve or by translating it in frequency. 

Involves the parameter p while the second Involves z^. 

Which should change? 

To distinguish between them, the program applies a 
test. If changing p forces a limiting curve to cross the 
measured curve, p will not be- changed . Then, z^ is given 
an opportunity. If can improve the fit without crossing, 
it will be changed. 

/ 

/ 
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PROGEAM 


Before beginning, the program must have the 
emissiyities , their frequencies and any fixed parameters. 

For each of the remaining parameters, a range must be 
established. The limits are grouped into two sets which 
generate the upper and lower limiting curves. These sets 
of limits must satisfy two requirements. First, these 
curves may never cross the measurements; the upper curve 
must be higher than the measured emissivity at any frequency 
Second, holding all other parameters constant, a limit 
from the upper curve, must generate a curve higher 

than the one generated by its lower curve counterpart, C 2 [j] 
as a parameter varies from the upper curve limit to the 
lower one, the curve becomes lower. 



figure 4 

■ .The model operates on both limits of each parameter 
to decide which better fits the measurements. Meanwhile, 
the other parameters are held constant midway between 


their limits. To -determine which limit fits better, the 
difference between the predicted emissivity, EMI, and 
the measured emissivity, em, is squared and summed over 
the five .frequencies . 

5 

VARI = ^ ( EMl[i] - em[i] 

i=l 

The limit with the greater variance is changed. In 
matching the measurements, the program minimizes the squared 
error. 

How much the limit should change is uncertain. Since 
emissivity depends nonlinearly upon each parameter, a 

e 

parameter range cannot be decreased by any set fraction. 

If, for example, the lower limit of the ice thickness, d, 
were moved halfway across the range, the actual value of 
d would be left outside rhe new range. Each change must 
be checked to see that the new limit fias not overshot the 
true value. 



73 


k limiting curve should be able to approach the mea- 
sured curve without ever crossing it. If the curves do 
cross, it would be because a parameter had bypassed its 
true value. This is the assumption behind a test to in- 
sure that the actual parameter value is still included 
in the new range. Before the old limit is changed, the 
proposed limit is tested to see that the curve it generates 
is everywhere above or below the measurements. If so, 
the sum of the absolute differences between the predicted 
and measured emissivities will equal the absolute sum 
of the differences. 

5 .5 

t 

i«l i=l ' . 

If this crossing condition is met, the new limit is adopted 

If not, the proposed limit will be moved back toward its 

♦ 

old value. 

The tentiV.i’iva limits change hy incrementing the 
index L, If L were one, and the limit did not satisfy 
the crossing condition, L would become two, and the new 
limit would be tested. 

Often it is not desirable to change the parameter 
range too drastically. The new limiting curve might be 
moved so close to the measured curve that subsequent para- 
meters' could not change without causing the curves to cross 


^ ^ ( EMI [i] - em[il ) 
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The initial value of L is called G and- must be 

specified each tim.e the program is run. When G is high, 

the limits- change gradually, and the. last parameter has 

much the same opportunity to converge to its true value 

tl tl 

as the first one has. is the j element in the 

first set of parameter limits. C2[j3 is its counterpart 
in the second set. 



figure 6 

The same way an area may be better filled using small 
squares, a curve may be better fit by approaching the. 
parameter values in small steps. 
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if k=l 


X=(C3^+G2)/2 


if k=2 


EMljiu] 2 s. 

vari=^ (eitt[m] -EMl[ni] ) ^^+var i > 

m : m+1 

J if m^5 j 

m=l ~ ^ 

if k=2 •*- — k=k+l 

© -«s : 

if- vari2>v 

then XjC^ ...i vari= ^ then XjC 2 


(Cj [j] •L+G2[3l) 


sumab=0 


if m<5 


if- vari2>varil 
then XjC„ 


X[j] = 


(C^[i3+L*G2[jl) 



^ EMXl[m] 

. sumn-eitt[m] -EMXl[m3 +suinn 
sumab = eia[m] -EMXl[m]+sum.ab y 

vari= (em[m] -EMXl[m] ) ^+vari 

m:m+l jT 

if m^5 ^ 

if sumabs^Isumnl 

if varil>vari2 

then C^[j]=X[j] else C2[j]=X[j] 


J:j+l-= 

^ n: n+1 

if n<^Y J 

C=CC^+C-)/2 

1 j2 

Display (C^,C 2 ,G) 


if 


FLOW GfiART 


PERFORMANCE 


I£ this program Is to be used with confidence^ the 
Inverted values It returns must be unique. No matter 
what the Initial limits are or- In what order the parameters 
change, the final result should be the same'. 

To examine the effects of the initial limits on the 
results, three sets of limits are inverted for the same 
set of emissivities . The first set. A, is. symmetrical 
about the true parameter values. The second set. A’, is 
asymmetrical. The third. A*', is identical to A' except 
for the range -of d, ice thickness. In A", the true value 
of the ice thickness is left outside the initial'* limits . 

For the set of parameters 


el=3.0 e"=0.05 p=l.-0 z^=0.i .. d=10, 

the Composite Model predicts these emissivities 


0.483 

0.537 

0.607 

0.705 

0.853 


0 1.4 GHz. 

0 2.7 GHz 

0 5.0 GHz 

0 10.7 GHz 

0 37.5 GHz, 


Initial Limits 


Final 

Limits* 

midpoint 

A 2.5 < 

e * 

< 

3.5 

3.21 

< 

e ’ 

< 

3.33 

3.27 

0.1 > 

e" 

> 

0.01 

0.06 

> 

e" 

> 

0.05 

0.054 

0 < 

P 

< 

2.0 

1.02 

< 

P 

< 

1.27 

1.15 

0.05< 


< 

0.15 

0.13 

< 


< 

0.14 

0.14 

20 . > 

d° 

> 

1 

12. 

> 

d° 

> 

•9.6 

10.7 

A* 2.0 < 

e « 

< 

3.5 

3.05 

< 

e ' 

< 

.3.23 

3.14 

0.1 > 

£»» 


0.01 

0.06 

> 

e" 

> 

0.05 

0.056 

0.5 < 

p 

< 

2.0 

1.09 

< 

P 

< 

1.27 

1.18 

0.09< 
15 > 


< 

> 

0.15 

-1 

0.13 

11 

< 

> 


< 

> 

0.14 

9.2 

0.14 

10 


*Twenty iterations were performed with G ^I^aI to nine 
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REPRODUCmiLITY OF THE 
ORIGINAL PAGE IS POOR 


Initial Limits 


Final Limits 


midpoint 


A" 2.0 < < 3.5 2.78 < e' < 2.96 2.87 

0.1 5 E» s 0.04 0.07 > e" s- 0.06 0,068 

0.5 < p < 2.0 1.08 < p < 1.27 1.18 

0.09<-2. < 0.15 0.12 < 2- < 0.12 0.12 

9<d^<l 8.<d^<7 ,7.6 



log }J 
(GHz) 


The final limits of the first two sets. A* and A", are 
almost identical. They differ in the value of e ' , but 
the difference is small compared to the wide range across 
which^'S' was originally allowed to vary. In this example, 
the midpoint values of A* are as good as or better than 
those of A. It seems there is no reason to believe that 
the initial limits should be symmetrical about the actual 
parameter values. 

Three of the five parameters have been very accurately 
determined. The midpoints for e ", p and d are in error by 
only ten percent of their initial ranges. Values for e ^ 
and Zq, however , ..are much less accurate. For A, the midpoint 
of ^ * is too high by more than twenty-five percent of its 


L. 
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initial range, is off by forty percent of its original 

range. For the error in the midpoint of s* is only 

ten percent of the initial range, but is too high by 

almost seventy percent of its initial range. 

This is the problem referred to earlier. During the 

inversion, a parameter is assumed equal to the average of 

its current limits. If this average is very different 

from the actual parameter value, the minimum error may 

occur when e’ is far from its true value. The program 

may be tricked into changing the better of the two limits. 

If the crossing condition then fails to detect the mistake, 

the new limit may leave the actual value of e* outside the 

range. For this curve, it is not surprising that e' and 

Zq overshoot their true values. and tend to their 

upper limits in order to match the very low emissivity at 

low frequency, and there is no point a limiting curve is 

1 

likely to cross as these two increase. 

The final limits for the set A" are quite different 
because of the incorrect range for d. Among the parameters, 
only E' has a midpoint value correct to ten percent of its 
initial range. Low as it was to start, with, the upper limit 
of d still decreased. .There is no indication that this 
limit was improperly set. When establishing initial limits, 
it is safer to make the ranges broad. - , 

So, for any set of limits that includes the actual 
value of each parameter, the results are the same. The 
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inversion should also be independent of the order in which 
the J)arameters are changed. The inverted results should 
be the same whether e* is considered first or last. Starting 
with the set of initial limits. A’, the same set of emis- 
sivities is inverted with the parameters considered in two 
new orders. Again, twenty iterations are performed with 
G equal to nine. 

The results for the three orders are almost identical. 
Again e* and are too high. Again e", p and d are 
very close to their true values. The only difference 
between the results is in e ' . There the difference between 
the midpoints is still small compared to the initial range. 
When G is equal to a large number like nine., the parameter 
order has little effect on the results. 


order 


Final Limits 


inidpo in t 


E", p, Zf^, d, e* 3.07 < < 3,25 ' 3.16 

0.06 > E" > 0.05 0.055 

1.09 < p < 1.27 1.18 

0.13 < z. < 0.14 0,14 

all.. > d > 9- CIO. 

d, z , p, £*•, £• 3.25 < £* < 3.50 3.37 

0.07 > £*' > 0,04 0.54 

1.04 < p < 1.22 1.13 

0.14 < z < 0.15. 0.14 

15, > > 6,7 11 

£*, E", p, z , d 3l05 < £' < 3.23 3,14 

0,06 > £" > 0,05 0.056 

1.09 < p < 1,27 - 1,18 

0.13 < Zf. < 0.14 0.14 

11 > d^ > 9.2 10 
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In this example, the inversion program determines 
e”, p and d very accurately, e' 'less accurately and fails 

I 

to determine For another curve, the program determines 

p and 2 q very accurately and e” and d less accurately < 

It seems that how well a particular parameter is determined 
depends upon the curve. In the first curve, there are very 
low values of emissivity at 1,4 and 2.7 gigahertz which 
force £* and Zq to be large. Because e’ and can increase 
without causing the limiting curves to cross any measurement, 
the crossing condition does not restrain them from over- 
shooting their true values. In the second curve, there 
is a distinct drop in the emissivity at 19.4 gigahertz 
that pinpoints Zq. There is- nothing, however, to highlight 
d and e". Their inverted values are less accurate. 


p = 0. 

75 

1 

z 

o 

II 

o 

0.776 


1.4 

GHz 

0.863 

0 

2.7 

GHz 

0.893 

@ 

5.0 

GHz 

0.873 

@ 

10.7 

GHz 

0.858 

@ 

19.4 

GHz 

0.870 


37.5 

GHz 


Initial 

Limits 

Final 

Limits 

midpoint 

3.5 5^ e* 

> 2.0 

2.72 > 

> 2.54 

2.63 

0.01< e" 

< 0.1 

0.07 < 

£" < 0,08 

0.075 

2.0 > p 

> 0.1 

0.98 > 

p > 0.75 

0.86 

0.15> z 

> 0.05 

0.073> 

z^ > 0.061 

0.067 

‘.1 :< d° 

< 50 

42 < 

° < 48 

45 
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INVERSION of REAL DATA' 


To test the model itself, emissivity data can be 
inverted to determine the physical properties of a site. 
Then the deduced properties can be checked with what is 
know about the site. 

T, Schmugge published an article including both mea- 

3 

sured emissivities and descriptions of three frozen lakes. 

Masato Nagase inverted this data using the Composite Model 

and correlated the parameter values he found with the 

4 

properties of each site. The inversion process was long 
and tedious. Using this program, the inversion can be 
done in minutes.* The results agree closely with Nagase*s 
and generate emissivities that differ from the measurements 
by less than 0,04 on the emissivity scale of zero to one. 

In performing the inversion, the real and imaginary 
parts of the dielectric constants for the underlying water 
or moist soil are assumed to be known. The ice thicknesses 
for the thicker two lakes are also assumed. Emissivity 
approaches a constant value past a certain depth beyond 
which the program cannot distinguish different thicknesses. 

• Eor Bear Lake, the thinnest of the three, the data is 
inverted for ice thickness. The ice is specified to be 
betweei one and fifty centimeters thick. The program 

narrows this range to thirteen to twenty-four centimeters. 

« 

*Eor five frequencies and five parameters, the time required 
per iteration Is about forty seconds. 

^ % 

eeproducibility of the 
82 ORIGINAL PAGE IS POOR 


The actual thickness measured at the site was fifteen 
centimeters. 

The emissivity data is given for six frequencies 
between 1.4 and 37.5 gigahertz. The emissivity at 37,5 
gigahertz is measured separately for horizontal and vertical 
polarizations. For the sake of matching, the two points 
are reduced to one, their average, 

Schmugge is unclear as to how the emissivity is measured 
at- 19.4 gigahertz. At one point, he refers to .that channel 
as horizontally polarized and scanning, In another place, 
he implies that that radiometer views from nadir, directly 
downward. Because the horizontally polarized emission is 
less intense than the vertically polarized emission when 
measured at an angle, the very low emissivity at 19.4 
gigahertz suggests that this radiometer does scan. Since 
the model works only ’for observations made from nadir, the 

-t 

measurement at 19.4 gigahertz is ignored. 

The initial limits assumed here and earlier are 

chosen to span the 'range of values found in ice and snow. 

The dielectric proper ties ‘ are reviewed by S. Evans in the 

■Journal of Glaciology.^ The parameters p and z^ are 

scattering parameters introduced by GurvichZ^ The choice 
* 

of range for them is based on Nagase^s experience which is 

7 

documented in his thesis. 
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frequency (GHz) 


aramet er 

midpoint 

Nagase * s- 

EMI 

em 


3.16 

results 

3.2 

0.497 

0.459 

i." 

0.034 

0,04-0.048 

0 .556 

0.582 

P 

1,39 

1.1-1. 6 

0.624 

0.608 

^0 

0,083 

0.08-0.12 

0.702 

0.720 

d=*15 

ewp=80 

18 


0.813 

0.800 


ewpp— 80 


Tae circles represent: the emis sivities predicted 
using the midpoint parameter values . 

figure 9 




parameter 


midpoint 


Nagase ^ s* 
result s 


frequency (GHz) 


EMI em 


1 

2,7 

2. 9-3. 2 

0.833 

■0.819 

tt 

0.082 

0.095-0.11 , 

0.905 

0.930 

p 

0.39 

0.35-0.45 

0.922 

0.893 

^0 

0.058 

0.05-0.065 

0.90. 

0.923 

d=50 



0.901 

0.899 

ewp=75 




■ 

ewpp=10 






figure 10 






frequency (GHz) 

arameter 

midpoint 

Nagase ' s' 

EMI em 




results 



e» 

3.0 

2.2 

0.8 08 

0.798 

e" 

0.049 

0.04 

0.887 

0 .883 

P 

0.82 

1,. 3 4 -1.3 8 

0.910 

0.925 

^0 

0.026 

0.037-0.04 

0.897 

0.914 

d=80 



0.847 

0.828 

ewp=7 0 


• 



ewpp=40 





The initial 
2.0<e’<3, 5, 

limits for 
0.1^e'*5«0.01 

these three inversions 
, 0<p<2.0, 0.0Kz-<0.1. 

were 

4 


^•dgur.e 11 
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SUMMARY . 

This general inversion scheme uses an iterative 
approach to tighten initially wide ranges about the actual 
values of a set of parameters. Limiting curves above and 
below the measured points are generated by two sets of 
parameter limits. These sets are made to converge subject 
to the condition that their limiting curves never cross 
the measurements. This criterion is intended to insure 
that the limits do not bypass the true parameter values. 
Because the parameters are considered sequentially, their 
order may affect the final results. To minimize this, 
the limits can be made to change gradually allowing each 
parameter equal' opportunity to converge to its proper 
value. 

The success of the program varies between sets of 
data. Some are better inverted than others. The quality 
of the inversion can always be checked by substituting the 
inverted values into the model and comparing the prediction 
withcthe -.measurement s . The result, however good or bad, 
are consistent for different initial limits and parameter 
orders. In the final examples, the program demonstrates 
its utility by successfully repeating three difficult 
inversions of remote sensing data done by M. Nagase. 
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(Cl) -'set eny fixer! pareneters . * $ 

(C?.) ev.'p i ev.'pp s SC$ 

(C3) v<:[C.Sf5e9, 1 . 715 elO^ 3.UCelO, 6w2Se10, ?,356ell]$ 


(C^i) "Set the data to be natchec'. , , ^ 

Ti is set cf enl ssI vi tl es v/as . generated by the 
Composite rode] for eip=2.f, eipp=C.09, P“1.5, 7-0=C.v..o'5, d-O.lS 
Wete that the unit of I enfth used Is the meter.' 5 


(C5) enii[0.602, 0.711, C.302, 0.827, 0.732]$ 


(C6) “Hcreis a key to I nterpretl ng the pararnter array, 
Cci p, ei pp, P, zO, d] 

Set the I r.i tl al 1 Imi ts, ■■$ 


(C7) ci:[2.0, 0.1, C, 0.0001, C.2]$ 
(CO) c2:[5.5, O.Cl, 2.0, 0,001, 0.Gl3$ 
(C9) "Catch the procram, ”$ 

(CIC) batch (zchu, 29, csk,konc) j 




(Cll) fit(c.y}:='jlcckc:<ecpflcat:true. 


ZERC, 
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L:G, 
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C:(CUC:)/2, , ^ 
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Kiip: EV(;;EALPART(0::..T'CM i*c[?] )) , :iVAL, 

E^I.:/CPA: 7(SC. T(^[1.]-K:i^CL2j)), EVAL, lAJAEA 

RGi: EV(3UAi(((?‘iIF-1 ) 

((MIr-!-l)^v-rH- EVAL, f.::AER), 

RI 2 : EVC3^AT( ( *2-i-(\IPr-^i.'P^*"2)/ 

( ) **2-!-(riIFr +HTrF ) *“2 ) ) , EV/l, t^L; :iR) , 
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S: 

E 
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(1 
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*9 
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(1-ER 


FF 

tvr 

L 



]/ 

3.F0 
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Ei-;i 

• 

. EV(1- 

.\-S, FV 

;.L 

, 
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. ' 1 

r»^ > 
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D I f 

p 

* 

r.\' i z.l 

. J -i 1 

t 


1, 

Mo: 

trr> 

V, 4 ^ 

), 



SCO 

p 

p 

EV(3I 

- * vV 7 

V 

L, 

Si 

i.ER) 

9 

), 



VAf> 

1 

: Et(V 

Ais I +SC.' 

9 

EVA 

L, 

MUI. 

r- “J 

111 \ 




:>i, 


IF i.<=-5 T..EN GCCb) ELSE r:i , 

IF R"! T-.EW VAAIi:7/,;I ELSE VARI2:VARI, 


1 

i\ ^ • I 9 


IF R=2 TPEli GO (A)., 

1 , 


I- v,*r.ii>VA,,i2 T.iTt^ i 3 'ir(r,i,!:]) else :;:^r-.’’EN:{Er.[]), 
IF \Yr,Il>v'F,.I2 T.iEW . j; : (L- F 1 [F2+Er[ J j)/(L-:-1 ) 
cirr rr.n : (T 1 r •- )/CL-^l ) 


.C, 


F, 


( *2 ; j , ]\ 

RXi?: [£V(oCi;;T ( ( ) **/ + (rijRRX->iTPP)**; )/ 
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SUr;W:Si:r.A2;:C, 

ax: !!y(rxi 2 *EXPC- 2 *t 4 :-p:*x[p]*v;[;!]/ 3 ,E 8 ), 

Gx: EV(r^^MJrX*x[43*;;r;;i/3.Er, mumea) , 

RX : EV ( (RX0 1 **2+AX )/{] + ( r,XO 1 *AX ) ^^*2 ) , MU; 'ER) , 

SX : E V ( ( ;; [ 3 3 *cx /^o * ( 1 ~ . ■ ; C 1 ** 2 ) * 

(l-ZXP(-A*'’lPPX*X[t J^l-:CMj/3.E3))/(i+Q;<:v^:2), HUXER), 

EHXi: EV(i-RX-sx, 

DIFX: EV(E:;X1-Er!;x], flU'ER), 

SUM ; EV(su:; +bjFX, MuxrR), 

SUMAR: EVCSUMA3-s-A35(DIF: ), HUMER), 
m : r.+i, 

IF K<=5 T.'iEti GCO) r.SE MM, 

IF SUM;M-'=A3S(SU/:ti) THEI^ G'^(E), 

IF L>1Q T.iE>l G0(F), 

L : Lt- 1 , 

GO(C), 

IF VARn>VARI2 THE C1M ELSE C2:X, 

REMYALUt ( VAR I , VAR I 1 , V/ R ! 2 ) , 

L:G- 
J : 1 , 
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1 , 
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G0(2ERl) 
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DISPLAYCCl, 
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3ATCH DO E 
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Abstract 


Variation of soil characteristics with water content/ 
soil type and its electrical properties are investigated by 
passive microwave remote sensing, A theoretical model is 
developed to distinguish the parameters or assess the 
possibility of distinguishing these parameters from wave- 
lengths ranging from X = 100 cm to X = 1 cm-. 

Radiometer measurements are freguancy, angle and polar- 
ization dependent. The theoretical model consists of an N 
layer medium on which various types of moisture profiles are 
mapped simulating surface v/ater distributions in desert areas , 
The results are expressed in brightness temperature and 
are compared with available data from the Nimbus 5 satellite , 
The validity and limitations of the model are examined as 
well. It was foimd that moisture content in soils with sand- 
like composition was not as easily detectable as the results 
of controlled experiments showed. 
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Introduction 


The field of passive remote £.cnsing is at the same 
time very old and very yomg* Very young if we remember 
that the first radiometer was acknowledged by R, H , Dicke 
in 1946 in a science revue entitled "The measurement of 
thermal* radiation at microwave' frequencies" in v/hich 
Dicke explains that the behavior of a microwave switching 
circuit that compares the thermal radiation power emitted 
from an object with that emitted from an internal reference 
source. -And yet, a field veiry old if one takes notion of 
all the variety of applications, not all fully exhausted, 
that were derived from it. Without involving ourselves 
with the whole list we may just mention a few applications 
of remote sensing such as measurements of ocean surface 
dynamics, ocean salinity, atmospheric temperature; detec- 
tion of sea ice, storm cells, agricultural crops and soil 
moisture content. We shall be interested with the latter 
feature . 

Whereas meteorologists are interested in tracking the 
moisture content of soil over large areas to learn more 
about mass and energy exchange at the air-soil interface, 
and hydrologists in the same behavior to predict runoff, 
our more modest task is to try to determine its presence 
in dry areas of mostly sendlike composition. 
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standing between the almost rainfree Sahara desert 
and the rain favored regions of the tropical zone, the Sahel 
region meets the definition of a transition area, and as 
such, displays many of the characteristics we intend to 
study. 

It seems thus far that we are nowhere limited in the 
analysis of thermal radiation except within the frequency 
spectrum. Dicke makes a point of confining himself to mi- 
crowave frequencies, and yet remote sensing at visible and 
infrared wavelengths has contributed useful data. However, 
it was noticed that the surface of the earth could not be 
observed when obscured by cloud cover, nor could it be ob- 
served during nightime, except perhaps with “cryogenically 
cooled far-infrared sensors". Hence the "operational utility" 
of these sensors is limited since most of the earth’s surface 
is covered by clouds or in darkness. In contrast cloud cover 
is much more transparent at microwave frequencies which make 
microwave sensors more insensitive to v/eather or sunlight. 

The only possible unattractive feature of microwave 
sensing is the poor spatial resolution associated with it. 

The resolution of optical and infrared sensors x^ill generally 
be better then that of microwave sensors since high resolu- 
tion can be achieved only with antennas large relative to 
wavelength. This is not a real problem, however, if the 
regions of interest are widespread enough such as ours. 
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But more determining factors are involved in our anal- 
ysis (as we shall see later) such as the sensitivity of the 
microwave radiometer, which as given by R. H. Dicke is: 

AT = 

Tjj == absolute system noise temperature 
B - RP noise bandwidth 
T = postdetection integration time. 

AT is the variation in brightness temperature ,v7he re 
brightness temperature is related to the surface temperature 

by T„ = e T„ (e - emissivity at -the surface) . 

s s 

2U: even more important factor is the drainage dynamics 
of water in sandlike soils . To understand better the phe- 
nomena we will resort to concepts in soil mechanics. 

Among the various techniques used to determine the 
brightness temperature or correspondingly the emitted power 
from the earth, radiative transfer equations are often the 
ones referred to. However, great care must be taken in 
applying them, since they are applicable in certain cases 
only. On the other hand through the Fresnel reflection 
coefficient approach we find more flexibility and dynamic 
range. We shall nonetheless investigate the different pos- 
sible approaches and choose the most appropriate for our 
particular analysis. 
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The Dielectric Constant 


Since the properties of soils are our major concern, 
we are thereby inclined to look upon the characteristics of 
their dielectric constants which depend on frequency, water 
content, and temperature . 

Frequency Dependence 

The frequency dependence of the dielectric materials 
is linked to the way they react to an electric field. The 
charge carriers within them are capable of being displaced 
by the latter and give rise to polarization. There are many 

n 1 

kinds of polarizations'" •' such as electron polarization a , 
atomic polarization a , orientation or (dipole) polariza- 

a 

tion a^, space charge or inter facial polarization a, and 
resonance polarization a^. We limit ourselves to the effect 
of the electronic displacement and orientation polarization. 
The Clausus-Mosotti equations expresses the relationship 
between the dielectric constant and the polarization as 
follows 1 



e + 2 3 


where n is the number of molecules’, and a = a + a ^ . 

e a. 

£ is the dielectric constant. 


page is po 
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When applied to moist soils the dielectric constant 
is allowed to be complex: e = e *' + i e". The imaginary 

part of the dielectric constant governs the dissipative 
feature of the soils. In Figure 1, tan 6 is the loss tan- 
gent tan 6 = sVe' arid t is the relaxation time (in the 
case of electron displacement alone it can reach the order 
of T Ip”^^ s) , 

The dependence of the dielectric constant (soils) on 
frequency is related to the dielectric relaxation spectrum 
of the water contained within it. 

When comparing the relaxation of water in bulk and 

(21 

moist soil. Figure 2 and Figure 3 , v/e observe similar 

distributions. The only difference is that "the maximum 
dielectric loss of relaxation in soils is displaced to lower 
frequency and the relaxation occurs over a narrower frequency 
band in soils than in bulk water" as noted by P. Hoesktin 
and A, Delaney in a paper treating with the dielectric pro- 
perties of soils 

An analytical expression was foxmd by Cole and Cole 
(1941) for the relaxation of soil water which is the followin 


“ e. 


e (w) = E 


CO 


1 + (iwx) 




a 


the subscript « and s refer to values of e. at frequen- 
cies such that WT >> 1 and 'wt << 1, t is the relaxation 
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time® a is a parameter showing the distribution of relax- 
ation times . 

When mixing water and soils without previously knowing 
the resulting value of ur it is easier to assume tlie sep- 
arate values of a (water) and a_{soil) and obtain the 

w s 

dielectric constant of the mixture by using 




e, 


ws 


- s 


s<» 


I + (iWT^)^ “ “w 


+ s ' + 

s>» 


ws 


- s 


S“ 


1 + (iwT^)^ “s 
s 


where t and t are the two different relaxation times, 
w s . 

a^, Cg represent the respective dispersion factor, (Actu- 
ally for water - !•) This equation suffers the incon- 
venience of not indicating the proportions in which the dry 
soil and water might be mixed. 

Water Content Dependence 

In that respect, mixture formulas are more informative 
and give us a better gauge on. the effect of water content 
on the dielectric constant.. It is useful, however, before 
reviewing some of these formulas to give ourselves some de- 
finitions regarding the terminology that will be used 
later. For instance, what is usually meant by moisture con- 
tent is the ratio of the weight of the water to the weight 


of the dry soil: water saturation is the ratio of the 

volume of water in the pores to the volume of the pores 
(if we assume sandlike soils then each grain of sand can 
be considered as having a spherical geometry* Within a 
configuration of four adjacent spheres the spacing in the 
middle can be assumed to be a "pore” whose size depends 
only upon the radius of the spheres) * Lastly, volume frac- 
tion of v^ater is just the ratio of the volume of water filling 
the pores to the total volume of the material. 

It is hard to find a single mixture formula that one 
can apply in all circumstances. Each one of them is only 
valid for either a particular frequency range, or mixtures 
having more than two components where the difference betv;een 
the dielectric constant may vary a great deal . 

Tarkhov's formula gives an average value of e for 
water bearing materials. 


E - 1 
e + 2 


e, - 1 e- - 1 e - 1 

^ 0 + _£ e + 9 

e- + 2-^ e^ + 2*^ e+2^ 

^ ^ T1 


where 0,, 6„f ... are the volume fractions of the com- 

1 A n 

ponents *** *’* ^ ®n ~ * ^1' ^2' '** ^n 

respective dielectric constants. For large differences 
between dielectric constants v;e have Lichtehicker *s formula 


log E = 0^ log + 02 log E 2 . 

For mixtures with only two components (soil, water) at 

. I 

microwave frequencies, we have Pierce’s formula 


E 



f^(l - F) 
1 - 




E„ ~ soil dielectric constant 
s 

~ water dielectric constant 
. = volume fraction of soil 

F = adjustable parameter 


But of all the formulas the one that showed the best ex- 
[ 5 1 

perimental matching was Wiener's formula 


£ 



S 


U + £ 

w 


U + f 
w 



(W) 


f„ = volume fraction of the two components f + f = 1 
s s w 

~ volume fraction of the two components f^ + f^ “ 1 

£^ = dielectric ■ • 

E — dielectric 
w 
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where F is an adjustable parameter 


u 


e + p 
w 


Temperature Dependence 

The third parameter, which affects the dielectric con- 
stant is temperature.' No correlation between the tempera- 
ture behavior of water in bulk and that in moist soils can 
be made since in the first case we have de V<^t < 0 and in 
the second de '/dt >0. 

However, a more quantitative relationship between the 

dielectric constant and temperature resides in the Debye 

« 

F6 1 

formula ^ , 


E 1 M _ 4 „ , 

^ • — = — Nit 1 a, 

£ + : 2 p 3 


+ 


u 

3kT 


M = molecular weight of the material 
p = density 
N = Avogadro’s number 

a = polarisability caused by electron displacement 
y = average dipole moment 
T = temperature in Kelvin. 


\ 
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Thermal Emission and Radiation 


D. H. Staslin, reports the physics of thermal emis'* 
sion as follows: *'A11 solids, liquids, gases aiid plasmas 

are coupled electromagnetically to their surroundings and 
hence absorb and emit electromagnetic radiation over a spec- 
tral region tanging from radio frequencies to optical 
frequencies and beyond. The emission results primarily 
from the thermal motion of the charged particles which com- 
prise the material, and depends upon the temperature of the 
material and upon the degree of coupling between the material 
and tlie electromagnetic field" . 


Blackbody Radiation • 

If we consider the case of perfect coupling to the 
electromagnetic field such as a perfect absorber blackbody, 
then the radiation intensity is given by Planck's radiation 
law: 

2 hf ^ -2 “1 -1 

X = watts m H ster (1) 

C' (ehf/M . 1, 


where 

h = Planck's constant 
f == frequency 


no 


k = Boltzmann ' s cons tant 


T « temperature. 

11 

Between the radio region and frequencies up to 10 Hz 
and for temperature greater than 5K (which is our case) , 
equation (1) can be reduced to the Rayleigh-Jeans Law: 


2 kT -2 —1 -1 

I ^ watts m Hz ster 


X = wavelength in meter. 


Prom this linear relationship between tlie intensity and 
temperature we get a similar relationship between the power 
P intercepted by an antenna (microwave sensor) and the 
temperature of a blacJtbody. We then have 

-1 

P = kT watts Hz 


or 


Tg = P/kB 

where B is the bandwidth of the frequencies considered. 
T_ is called brightness temperature , 


Emission by a Surface 


We can associate a brightness temperature to a surface 
(earth, air boundaries) . "If the intervening atmosphere 
between those surfaces and the radiometer is lossy then the 
radiation emitting from them will be somewhat attenuated 
and the lossy atmosphere itself will emit thermal radiation, 
A radiometer would thus yield an apparent absolute tempera- 

r 8 1 

ture Tg given by the radiative transfer equation" '■ , 

Tg = t(e + r Tg) + (1 - t)Tg (2) 

t = atmospheric transmissivity (0 t <_ 1) 

e = emissivity of object ® 5. 

= absolute temperature of the surface 

r = power reflectivity of surface £ I* 5. 

Tg = equivalent absolute radiation temperature 
incident on surface 

Tg = absolute temperature of lossy atmosphere between 
surface and radiometer. 

One of the interesting parameters of equation (2) is 

the emissivity e. By applying the principle of detailed 
[91 

balance we get that; (or reciprocity) 
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E “ 1 “ r 

r = I eP 

where if we consider a half-space such as in Figure 4 with 
an incident electromagnetic wave, R is the reflection 
coefficient and r is the reflectivity. R can be expressed 






Different Approaches to Surface Sensing by 
rtadiative Transfer Equations 

As mentioned earlier, the brightness temperature is due 
to emission of radiation of earth surfaces. However, the 
word surface is ambiguous . Is the surface the only part that 
contributes to the emission of the earth? Certainly not 
since the whole subsurface region adds to the total emission . 
To determine the depth at which the emission of radiation is 
negligible we refer to "reprocity" and eKamine again tlie case 
of an incident electromagnetic wave. The probing depth (or 
equivalently, the depth just mentioned above) would be the 
one at vjhich the electromagnetic wave would be totally at~ 
tenuated. , 

Ha If "Space Model 

It is through such an approach that K. F. Kunzi modeled 

the brightness temperature of snow and ice. Hopefully, 

Kunzi ’s model for snow and ice will enlighten us in soil 

moisture sensing. He uses the radiative transfer equation: 

i 

0£(x' , V )djc' \ 

’ dx . (3) 

Here T(x) represents the subsurface temperature, a(x, v) 
is the absorption coefficient (or attenuation coefficient) , 


( r- 

Tg(®k) = E(v)(y . T(x) a(x, v)e 
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X is the depth, v is frequency, and E is the surface 
emissivity. We neglect atmospheric losses. 

The form of a(x, v) can be found, if we consider 

the electromagnetic fields propagating in the x direction 

-i (ut - k x) 

wiiha dependence of the form e * where is 

the component of the propagation vector' in the z direction 
We shall assume normal incidence and take k„ - k. We fur- 
thermore have that k ,= , since snow and ice have also 

complex permittivities . k has then a real and imaginary 
part, which implies that the fields will decay exponentially 
with depth (Figure 5) . The attenuation coefficient a is 
such that 

a “ Im(k) , . • ■ 

Since the half-space considered is inhomogeneous, 
is a function of x, and since e is also a function of 
frequency, a(x, v) ~ Im(k) , 

The sursurface temperature on the other hand is approxi 
mated by 

T(x) = T^q + 6 (4) 

4 • 

= the temperature at the surface 
Y - the thermal decay constant 
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6 “ temperature deviation from 

\fhile neglecting scattering, the other assumption being 
made is that for the dielectric constant, e* is independent 
of temperature and frequency. It was also found that 
e' >> e”, hence that power absorption coefficient can be 
approximated' by: 


a 




2ttve» 

c/rr ‘ 


Integrating equation (3) with a held constant {vari- 
ation of with temperature and depth are neglible for 

ice and snov?) and using equation (4) we find 

* 

E ’^10 -f- ■ 

oi + Y 

Let us then try to apply Kunzi's approach, to the case 
of soil moisture determination. To simplify our task we 
consider the case of constant temperature throughout the 
half-space. Equation (3) can be transformed in the following 
way; 
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The problem here is to be able to determine the moisture 

content of the soil given the brightness temoerature T-, , 

s 

If we want to get a closed form relationship between 
Tg and the moisture distribution in the soil f(5c) we see 
from equation (5) that a(x) must be integrable for this 
to be possible. 

But what is the dependence of a on f (x) . Using 
Wiener's formula; 


e “ 


f„ U + e f 
w w s s 

f U + f 
w s 


V7here 


e + F 
U := , 

+ P 




identifying f^ = f (k) . 


e,. = dielectric constant of water 
w 


e = dielectric constant of soil . 
s 


We obtain, 


e = 


U f(x) + Sgd “ f(x)) 
f (x) U 4- (1 - f (x) ) 


e = 


• U “ eg)f(x) + Sg A f(x) + B 


1 + (U - l)f(x) 


C f(x) + D 
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A, Br C, D, are complex. Therefore: 


a(x) = ^ Ira . 

c V C f (x) + D 

Suppose j we are given T_, T„, w, e., e„. V'7e find 

that we cannot determine e uniquely since there are an in- 
finite number of ways to relate a "real*' value of the square 
root of a complex number . 

We could on the other hand postulate the form of f (x) 

and tentatively (a computer helping) try to approach the given 

With 

value of Tg, but in doing so we may loose touch of the phys- 
ical reality of the moisture distribution in the soil. 

The other deceiving feature of this technique is that 
tlie constraints imposed on Kunzi's model make its application 
to moisture sensing less accurate than what might have been 
expected. This is due mainly to the fact that the attenua- 
tion (moist soil) " was " significant compared to a wavelength, 
and because the real part of the dielectric constant changes 
rapidly . 

Generalized Half-Space Model 

To remedy this, we hope that a more general approach of 

ril 1 

radiative transfer such as stogryn 's would be more easily 



applicable to our case. The procedure Stogryn uses is still 
with a half“space model with vertical variations in permit- 
tivity and temperature. He finds a radiative transfer equa- 
tion by solving Maxwell's equations using a randomly fluctua 
ting current density with a known mean value , Since his 
development is a rather lengthly analytical discussion of 
Green's function solutions with boundary condition to be 
satisfied we shall go, directly to his results which stand as 

Tg = 1 r| 2 + [1 - [RjM f dx' T(x') f(x’) 

•/ ■»DO 


where 


o))3 / JSP J- \ exp l“2Im f ydx"] 

■IvCx'jJ 


Re {(e(0, w) - sin^0) ' } 


= normalized dielectric constant 
= Cw/c) ( {X| (ij) - sin^e)^^^ 

:= angle of incidence of electromagnetic wave, on the 
half-space (we shall take its value to be zero in 
our case) 

- can be considered as the surface temperature T^ 

= temperature profile in the half-space.- 


Stogryr^'s model is very appealing even though the 


REPRODUCIBILITY OP THE 


A T, _ JEXA 



expression of brightness temperature is long. However, 
because it is limited by the constraint of a slowly varying 

j 

dielectric constant (this time real and imaginary parts can 
vary) 


ds 

< < 

(1) 

— £ 

dx 


c 


this would limit in turn our analysis of moisture sensing 
to slowly varying distribution which may not always be the 
case . 
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Skin Depth Approach 


Rather than using the radiative transfer approach to 
the problem/ G. Poe and A. T. Edgerton relied on Fresnel 

ri2 1 

Reflection Coefficient for a half-space as in Figure 4. 
The model they use to ascertain the moisture content in soil 
is as follows: "the model originates from the hypothesis 

that it is the total mass of water per unit area lying be- 
tv/een the surface and one e-lectroraagnetic skin depth which 
determines the brightness temperature". 

The relation relating the skin depth to the dielectric 
constant of the half-space/ is: 


L ^ 1 

X 2ir Im 


where 6 is the skin depth, X the wavelength, k = 
with the familiar eKpression of the emissivity 

e - 1 - r ‘ (7} 

and Wiener’s mixture .formula (W) they derive a four step 
procedure to obtain the total mass of v/ater per unit area. 


EEPRX)DUCIBILI^^0I 

OEIGINAb PAGE 13 '■ 


123 


is the 


(i) Assuming that' the surface temperature T 

same throughout the medium, then from measured brightness 
temperature (corresponding to specified volume percent mois- 
ture content) and an effective emissivity e = 

is determined, 

(ii) From the established emissivity in part (i) and 
the corresponding volume percent of moisture content, a 
curve relating e and percent of water in the soil is drawn, 
(Use equations (7) and (W) .) 

(iii) Then a curve relating 6 to percent of water is 
determined with equations (6) and (W) . 

(iv) Finally, the total mass of water per unit area is 
obtained .by just multiplying 5 and percent of water. In 
the result of their analysis it is found that the brightness 
temperature behaves generally linearly with moisture con- 
tent. (Figure 6.) This is in accord with some experimental 

ri3i 

results obtained by T. Schmugge 
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Proposed Model and Approach 

We gave ourselves the task of determining the soil 
moisture content in the dry region of the Sahel, However, 
unlike Poe and Edgerton, the emphasis will be made more on 
determining a moisture profile rather than moisture in bulk . 
Trying to be more qualitative implies a more thorough know- 
ledge of the dynamics of water flow in sandy soils. The 
problem at hand consists in evaluating tlie shape of the water 
content profile which would simulate the penetration of v/ater 
after rainfall , 

Soil Mechanics Considerations 

The drainage mechanics in cohesionless soils has been 
fully explored by T. W. Lambe^^^^. Hence we shall appeal to 
his theory for the understanding of water flov; in sandlike 
soils. In a first step we refresh our memory about the gen- 
eral hydrodynamic equation of flow. The assumptions are, 
as in Figure 7, that flow. is vertical and laminar. It will 
be laminar if we assume that the soil is made of fine sand, 

We are within that limit if we consider a grain size dis- 
tribution such as in Figure 8. By using Darcy's Law we ob- 
tain the rate of flow r in a soil element such as in Figure 
6, V7e then determine the net rate of flow, Ar by consid- 
ering the rate of flov; both at the bottom and at the top of 
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the soil element (z direction) . We then equate Ar with 
the internal change of water within that element and obtain 
finally the general hydrodynamic equation for laminar flow 
in soils which is the following: 


k. 


r,2-u 3k„ ah 

o rx ^ 


3s‘ 


az‘ 


ak 


3z 


2-dz 

3z^ 


3G „ 3e 
e — - + G <— 


1 + e 


at 


at 


( 8 ) 


where 


k = permeability in the z direction 
2 

h - head 

e = V is the void ratio 
G V /V is percent saturation . 

If V is the total volume dx dy dz then, 

V = volume solids + volume voids - where 

V = volume water + volume air = V + V . 
v w a 


V =s dx dy dz 

^ 1 + e 


Ge 

V = dx dy dz 

^ ■ 1 + e 
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V = dx dy dz ■ 

® 1 +.e * 

♦ 

If in the above equation k is constant, e is constant, 

G is constant then we get 8^h/8z^ = 0 which is the fami- 
liar Laplace's equation for one dimensional steady flow. 

When considering flow in soils we can generally focus 
on four groups , they are as follows ; 

(i) Flow where G remains constant but e varies , 

(ii) Plow where e and G are both constant or 
steady flow, 

(iii) Flow where e remains constant and G varies. 

(iv) Plow where both e and G change. 

^Che two first groups have been thoroughly studied 
and represent no great interest here, Lambe focused on the 
third category which is no doubt, of greater interest to us. 

The important aspect as shown by Lambe is that flow in 
solids depends upon a total head h, which in turn is made 
up of several types of head. One type of head which con- 
tributes greatly to flow is the capillary head. By head we 
mean here a cannular path which is used by the water to flow 
within the soil, iflhen discussing movement of v;ater in soils, 
we must mention then the effect of capillary pull since there 

\ 
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is a tendency for the water to be held np against gravity, 

Lambe experimentally verified this phenoma. The curves 
that were obtained relating height and corresponding level 
of saturation agree with the following equation: 


-1/2 (H /x')^ 

S(x') = 1 - e ® . (9) 

where S(x’) is the saturation percentage from a height 
x' of total saturation. 


II = capillary head and is defined as 

TT - 2T 

gpcr 

where T is the surface tension, g is gravity, p is 
the density of water, o is the void radii which will be 
generally considered as Rayleigh distributed. 

V7hat we mean by Rayleigh distribution is a behavior of the 
probabilistic distribution of the void (pore) radii as in 
Figure (9) , If the probability distribution of void radii 
r, P{r) is such that. 


f r g-r /2b ^ ^ 0 

b 


P{r) 


0 otheirwise 
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(where b is a constant and varies with soil type) then 
the radius associated with the expected value of P (r) can 
be roughly taken as our a in the definition of , The 
expected value of P(r) is taken to be: 

P(r) P(r) dr - . 

J -0 ^2 


Fresnel Reflection Coefficients for a U-Layer Media 

Now that we have more insight on the soil mechanics 
aspect of our study we can give ourselves a model that will 
relate a brightness temperature to a moisture profile con- 
forming to the one just investigated. We shall consider a 

half-space medium, witli constant temperature T , with no 

s * 

scattering, non-magnetic which is inhomogeneous in permit- 
tivity if we consider the mapping of a moisture profile onto 
it as in Figure 10. We then partition this half space into 
N layers of uniform permittivity (calculated through Wiener's 
formula) . We have to then reduce the problem to finding the 
reflection coefficient of a N layer stratified media. We 
use Fresnel reflection coefficients as reviewed by J. A. Kong 

risi 

■* to find the total reflection R at the surface (air, 
soil boundary) . 
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f 


R = i 


- 1 / 

“ •** /— ' + ^^01 d^) 

^OX / ^01 


L- + ^2 


®12 / ®12 


lU 


+ EJ3 expI12k2^(d2 - a^)] pS . ' . y ^ ... 


^23 / ^23 


\ > ^ 


■" ^{n - l)n ^ ^nt ^ ^ 


> . 


( 10 ) 


where ^(n - i)n reflection between two adjacent 

partitions as we go deeper in the sbil* 

d„ is the width of each uniform layer? k is the 
propagation vector of the electromagnetic wave in the x 
direction, (Actually we may just call it k^ since we are 
considering normal incidence of the E and M waves on the 
medium.) 

To determine the emissivity we recall, 


E = 1 - |r|^. 


( 11 ) 


With Tg taken to be the temperature of the whole subsurface 
region, we then obtain the brightness temperature. 
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Tg(®k) == E • Tg,. 


( 12 ) 


Through equations (10) , (11) , and (12) , we can establish 
curves relating E (and ultimately T^) to frequency for 
various heights of the moisture profile (in doing so we test 
the probing sensitivity of microwave sensing) , In view of 
equation (10) ^ the numerical analysis involved here requires 
a computer program such as the one listed at the end of these 
pages , We must however comment on the meaning of some of the 
characteristic variables used, 

HS, is the capillary height which determines the 
maximum height of our profile. 

KS, is the value of the dielectric constant of the soil 
(KS = 2,55 + i 00.015 at frequencies 3 GHz, 

KS == 2,53 + i 0,009 at frequencies above 3 
(normalized to Sq) f 

FREQ is the value of the frequency in each considered 
frequency reinge we are interested in (in GHz) be- 
tween 0,3 GHz and 31.4 GHz. The frequencies con- 
sidered are successively: 0.3, 1.0, 2.0, 4.0, 10.0, 

22 .2, 31 .4 GHz . ^ 

P is the form number used in Wiener's formula. It is 
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such that for frequencies < 10 GHz F =* 32 and for 
frequencies > 10 GHz F = 16. 

, • 

kw is the dielectric constant of water. 

z is the spacing between each layer of uniform dielec- 
trie constant. 

EX, its imaginary part, is the decay or attenuation 
(of E and M wave) factor. 

CF is a proportionality factor which varies with HS 
to give saturation at the bottpm of the head and zero 
at the top . 

We will notice that saturation, in the profile we use 
is not 100% but 40%. This is due to the fact that in 
Wiener's formula the water mixture has to be in volume per 
cent. Therefore, if saturation by weight is considered to 
be 20%, and our type of sand has a specific gravity of 2,66, 
then the saturation level in volume percent is roughly 40%. 
Lastly, N, is the number of layers taken in the model. The 
value should be dependent upon the height of "transition zone* 
Ax as in Figure li where the change in the dielectric con- 
stant is the fastest. Another criterion for the choice of 
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N is related to the frequency. We assume that K must at 
least satisfy N = 10d/3X where d is the depth of the 
profile and X is the wavelength. Generally, for better 
sensitivity and accuracy (as to the form of the profile) we 
take ten to fifteen times the value of N obtain by the 
above criterion . 

There are three distinct ways in which the profile can 
vary: the height of the profile, the moisture percent at 

the surface, and the moisture percent at saturation. 

Computer Results 

We shall first consider variations with height of the 

moisture profile. The heights examined were successively: 

3^ = 0,7 cm, ^2 ~ 1.5 cm, 3^ - 5.9 cm, 11.0 cm. 

Beyond , the lowest frequency (our reference point since 
* 1 

it has the greatest penetration depth) gives invariant re- 
sults in emissivity (E .95 always) . The results in 
Table M show that the most sensitive frequencies are the 
■three lowest ones. However, the amplitude of the emissivity 
variation is only in the order of 0.04, (This detection 
process would thus require an emissivity sensitivity in the 
order of 0.01 at least.) 

k 

We now vary the moisture percent at the surface of the 
soil (hence simulating recent rain) , The variations are 
taken from 0% to 20% as in Tables B2 (a) and B2 (b) . This 
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Frequencies 


0.3 GHz 1 

GHz 2 GHz 

4 GHz 

10 GHz 27.2 GHz 

31.4 

(1) 

(2) (3) 

(4) 

(5) (6) 

(7) 

E (n) V7xth 

n = 1 ... 7 

is the 

corresponding emissivity 

z - 

0 .-7 cm 

1 . 5 cm 5 . 9 cm 

11 cm 

E(l) = 

0,88 

0.88 

0.89 

0.92 

E(2) = 

0,88 

0.89 

0.90 

0 .93 

B(3) = 

0.88 

0,91 

.0,93 

0,94 

E(4) = 

0.91 

0 .91 

0 .94 

0.94 

E(5) - 

0.94 

0.94 

0.94 

0.94 

E(6) = 

0.94 

0 ,94 

0,94 

0,94 

E(7) = 

0.94 

0.94 

0.94 

0.94 



Table A1 

s 



Effect of height variations of the moisture profile* 


GHz 
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time the whole spectrum of emissivities is affected and 
furthermore by the same amoxmt (AE = 0.01) . 

However, in one case (Table 2.) the geight is 
and the moisture variation is from 0 to 20% and in the 
other the height is with a moisture variation from 

0 to 10%, This shows better sensitivity at lower depth. 

Finally’/ we vary the moisture percent at saturation, 
(Actually in the model the saturation level is uniform 
after a certain depth which v/e called z^ n = l, 2, 3, 4. 

It is at that depth, -that we make the variations in moisture 
percent.) Varying the volume (what we called moisture) per- 
cent at saturation reflects a variation in the specific 
gravity of the sand that constitutes the soil. We took a 
saturation percent of 40% which corresponded to a specific 
gravity of 2.66 for the sand grain. For the purpose of 
discussion v;e considered the case of a specific gravity of 
4 which corresponds to a saturation volume percent of 50% . 
(However, as mentioned, the characteristic sand would be 
considered extremely fine. This may therefore alter the 
shape of the moisture profile. We are therefore considering 
a hypothetical situation.) The results are sho\-m in Tables 
3a and 3b. Chaiiges are noticeable only in the case of the 
two lower frequencies-. Lastly, we consider the case of si- 
multaneous variations -as in Tables 4a and 4b where we varied 
in both cases all the parameters discussed so far. The 
variations in emissivity were still bound to AE = 0.04 for 
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Moisture content 0% 


Moisture content 20% 


E(l) 

£S 

0,93 

E(l) 

= 0.92 

E(2) 

= 

0 .94 

E(2) 

= 0.93 

E(3) 


0,95 • 

B(3) 

= 0.94 

E(4) 

= 

0.95 

EC4) 

II 

o 

* 

E(5) 


0 .95 

E(5) 

= 0.94 

E(6) 

= 

0 .95 

E(6) 

= • 0 .94 

E(7) 


0.95 

B(7) 

= 0.94 


Table 2 (a) , 

Effect of surface moisture variations on emissivity at 

z = 11 cm. 
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Mo is fcure con ten t 

E(l) = 0.90 
E(2) = 0.90 
E(3) := 0,94 
E(4) = 0.95 
E(5) = 0.95 
E(6) = 0.95 
E(7) = 0.95 

Effect of surface 


0 % 


^5oisture content 10% 


E(l) = 0.89 
E(2) = 0.90 
E‘(3) ^ 0.93 
E(4) = 0.94 
E(5) = 0.94 
E(6) =0 .94 
E(7) = 0.94 

Table 2 (b) . 

moisture variations on emissivity at 
IS = 5 . 9 cm . 
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Moisture content 40% 


Moisture content 50% 


E(l) = 0.93 
E(2) = 0.94 
E(3) = 0.95 
E{4) = 0.95 
E(5) = 0.95 
E(6) = 0.95 
E{7) = 0.95 


E(l) = 0.92 
E(2) = 0.94 
E(3) = 0,95 
E(4) = 0.95 
E(5) = 0.95 
E{6) = 0.95 
E(7) = 0.95 


Table 3(a) . 

Effect of bottom saturation level on emissivity at 

z = 11 cm. 
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Moisture content 40% 


Moisture content 50% 


E(l) = 0»90 
E(2> = 0.93 
E(3) = 0.94 
E(4) := 0.95 
E(5) 0.95 

E(6) = 0.95 
E(7) 0.95 


E(l) = 0.88 
E(2) ^ 0.93 
E(3) = 0.94 
E{4) = 0.95 
E(S) « 0.95 
E(6) = 0.95 
E(7) = 0.95 


i 

Table 3 (b) . 


Effect of bottom saturation level on emissivity at 

z == 5 .9 cm. 
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2 = 0 , 7 cm 


2 « 11 cm 


Surface moisture content 10% 
Bottom saturation level ■ 40% 

E(l) = 0.88 
E(2) = 0.88 
E(3) = 0.88 
E(4) = 0.91 
E(5) “ 0.93 
E(6) = 0.94 
E{7) = 0.34 


Surface moisture content 0% 
Bottom saturation level 40% 

E(l) ^ 0.93 
E(2) = 0.94 
E(3) = 0.95 
E(4) = 0.95 
E(5) = 0.95 

» 

E(6) = 0.95 
E(7) - 0.95 


Table 4 (a) . 

Simultaneous variations of parameters. 


\ 
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z - 1»5 cm 


Siirface moisture content 20% 
Bottom saturation level 40% 

E(l) = 0.88 
E{2) ^ 0.88 
E(3) == 0.91 
E(4) = 0.94 
E(5) = 0.94 
E(6) = 0.94 
E(7) = 0.94 


Surface moisture content 0% 
Bottom saturation level 40% 

E(l) = 0.90 
E(2) = 0.93 
E(3) =0 .94 
B(4) = 0.95 
E(5) = 0.95 
E(6) ^ 0.95 
E(7) = 0.95 


Table 4 (b) , 


Simultaneous variations of parameters. 



most frequencies. 

Three major features can be deduced from this analysis. 
First the probing depth of the detection process is limited 
to no more than 20 cm. Second ’the emissivity variations 
are limited within a small range AE - 0.04. Third the 
emissivity is dependent on the surface moisture. 

Possible Future Investigation 

If the latter varies with time (after rainfall^ water 
will filter tlirough the soil gradually) then the emissivity 
will also be a function of time. However/ the longer the 
elapsed time between rainfall and the probing procedure the 
less likely the moisture content (or distribution) will be 
distinguishable. W, T. Lambe derived experimentally that 
the moisture profile (related to the capillary effect) would 
vary in fact with a characteristic time T where 

^ _ k dh t 
n dG Ii"" 

where k is the permeability, h is a head, n is the 
number of heads, L is their average height, and t is 
time expressed in minutes. 

It w?as found that for a characteristic value of T 
(for which the surface moisture content goes from saturation 
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to 0% water) corresponded to an elapsed time t ~ 120 mn, 

fl71 

Also it was noted by T. Schmugge that the total 
moisture content within a particular surface layer does not 
change readily. He states "Results indicate that, although 
the surface layer changes rapidly (He means the distribution,), 
tlie total moisture in the top 15 cm does not: thus the 
average soil' in the top 15 cm decreases from 20% at 3 
days to 16% after 6 days, whereas the surface layer top 
half -"Centime ter has dried 'from 20% to 8,5% ", 

•We see that if we were interested in times in the order of 
2 hours instead of 72 hours then the decrease in the top 
15 cm would not be nearly as great as 20%. If, furthermore, 
we used satellites with orbital periods of 2 hours, then 
we might investigate a moisture pulse propagation in sandlike 
soils , 

It would be reasonable as a first approximation to 
assume a gaussian distribution of volume percent in the soil 
as in Figure 12, The purpose would then be to determine v^ 
the moisture propagation velocity given Tg(t) as a finction 
of time. The moisture profile would be of the form: 


M(x, t) 


-(x-vt)'-st 
e ^ e 


where = magnitude (maximum volume percent) 
Vp = propagation velocity 
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= smearing constant 
K » depth from the surface . 

The smearing constant s would in fact be a function of 

m 

porosity, temperature and humidity content. Of course no 
work has been done in that area and no formula has been given 
to determine 

Sj^(P, t, h) 0 


(Figure 12.) 


Having familiarized onrselves with the problem of 
moisture sensing it seems that we have a better chance of 
interpreting the available data related to tlie Sahel region 
in Africa. Our task would be somewhat pretentious if we 
v;ere to consider the whole of the Sahel, We have therefore 
confined our study to the following geographical sites: 

12“W, 4‘^W, 4®E, 12®E, 20°E of longitude and 17°N of 

latitude. These sites have the interesting feature of being 
located at tlie center of the Sahel belt. 

Data from Nimbus 5 

The satellite used to record data was tlie ' Nimbus 5 
which was launched in late 1972. Its microwave sensors were 
locked on the following resonant frequencies: 22.2, 31,4, 

53,6, 54.9 and 58.8 GHz, We shall' only be interested in 

the two lower frequencies since they lie within the water 
vapor and liquid water emission band respectively. The high- 

n 8 1 

er frequencies are much less sensitive to thera^ ^ (mainly 
used for probing the atmosphere) . 

The brightness temperature measurements were generally 
taken within a time period of seven to eight days at intervals 
of twenty days . This was motivated by the fact that the 
flight pattern of the satellite had to overlap on five specific 


locations. It should be noted that all the measurements 
were taken over a period of seven months (from December to 
June) . Some months have been omitted for lack of cataloged 
data. Furthermore, the readings in Data Set l of location 
(3) and (4) of the month of December bear the year 1973 
instead of 1972. This is due to the launching time of the 
satellite which created an offset in the space-time corres- 
pondence for locations (3) and (4) . Consequently, no 19 72 
data was given for them. To remedy this, readings were taken 
in December 1973, The locations (or sites) were chosen with 
at eight longitudinal degree interval to best eliminate any 
possible overlap of readings that could be caused by the 
200 km X 200 km resolution coverage of the satellite. 

Ground Truth Data 

However, to have a more complete status on. the region 
covered ground truth data is needed. Table (1) displays such 

[ 14 ] 

information . The stations from which the readings were 
taken are in accord with the locations of interest. 

latitude longitude 


Timbuctu 

16 '‘N 

5®E 

Bilma 

18“N 

IS^E 

Faya-Large au 

18®.N 

19 “E 

Nema 

16°N 

7°W 

Atar 

18®N 

13®W 
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yABIE-(l) 


I 


YEAR REGION 

JAN 

FEB 

MAR 

■ AFR 

MAY 

JUII 

JULY 

AUG 

SEP 

OCT 

NOV 

DEC 


\ 1968 TD 4 BUCTU 

295 

298 

301 

304 

307 

3O6 

504 

301 

304 

304 

301 

296 

(K) 

■ 

0.0 

0.0 

0.0 

0.0 

2.5 

18 

57 

80 

35 

2.5 

0.0 

0.0 

(mn?) 

1968 BIUyIA 

290 

293 

294 ‘ 

301 

305 

306 

306 

506 

304 

300 

295. 

290 

(K) 


0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

2.5 

12 

5 

2.5 

0.0 

0.0 

(irnr) 

1972 BILMA 

300 

303 

3O8 

511 

‘ 315 

316 

315 

515 

312 

511 

506 

300 

(K) 

. 

0.0 

0.0 

0,0 

0,0 

1 

1 

3 

10 

5 

2 

0.0 

0,0 

(iPip) 

1972 FAYA-LARGEAU 298 

303 

307 

312 

514 

515 

314 ■ 

515 

312 

310 

306 

301 

(K) 

1972 REMA 

505 

506 

510 

515 

515 

315 

511 

3O8 

310 

312 

309 

504 

00 


0.0 

0.0 

1 

5 

10 

50 

63 

• 111 

55 

15 

0.0 

0.0 

(eic) 

1972 ATAR 

500 

302 

307 

509 

311 

315 

515 

314 

315 

310 

306 

501 

(K) 


1 

1 

0.0 

0.0 

1 

8 

95 ' 

55 

55 

8 

10 

1 

(n?ir) 

REGION 

JAN 

FEB 

MAR 

APR 

im 

JUN 

JULY 

AUG 

SEP 

OCT 

NOV 

DEC 


■RIU/IA 

1 

1 

1 

1 

1 

1 

2 

2 

1 . 

<1 

1 

1 


FAYA-LARGEAU 

2 

2 

2 

2 

2 

1 

2 

■5 

2 

1 

2 

2 


WEt<IA 

5 

5 

2 

3 

4 

3 • 

4 

5 

4 

3 

5 

5 


AM 

3 

2 

2 

2 

3 

3 

3 

3 

3 

4 - 

4 

4 



( ffiGREE OF CLOUDINESS ) 
in oktas 


\ 
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The temperatures in Table (1) are expressed in Kelvin degrees 
since they will be used as ground surface temperatures. The as- 
sumption here is that the ground and the ambient air are in thermal 
equilibrium. The difference in the measurements at Bilma between 
1968 and 1972 (temperature and rainfall) shows clearly the trend 
towards a drought. We may specify that all measurements of tem- 
perature and rainfall are averages (unfortunately) per month. 
Cloudiness readings are inoktas^ where one okta correspond to a 
degree of cloudiness where l/8th of the sky is covered by clouds 
(mostly cumulus). Also, it has been recorded that evaporation was 
highest during the. pre-rainy season. But the mean relative humi- 
dity in July is approximately the same as that in December for the 
latitude considered (the level of humidity is invariant over the 
year) . 

As for the geological features of the Sahel, they are mostly 
uniform. We can give ourselves a sort of typical profile of la- 
terite soils if we examine the region of Sefa, Senegal ^ , 


Parent rock: 
0-12 cm; 

12 - 28 cm: 
28-55 cm: 

90 - 120 cm: 
120 < cm; 


Sandy-argillaceous sandstone 

Carboniferous surface debris, slight cohesion, porous 
Sandy, medium cohesion, porous 

Hardened, sandy argillaceous, strong to medium cohesion, 
porous 

S an dy -argil lace ous , unhardened concretions 
Transition to argillaceous sandstone. 


Of course, this profile will be siibject to variations 
as we go from the Western-most part of the Sahel to the 
Eastern-most part of it. However, generally we are likely 
to find sandlike characteristic in the soils. 

We have so far mentioned only the existence of a data 
set I. Actually there are. two data sets, I, and II. This 
is due to the confusion (which may be unexcusable) that arose 
when we tried to match our locations to the satellite's tra- 
jectory. As a result, in the first retrieval of data, three 
of the five locations gave emissivities which went as low as 
E « 0.51. This being more characteristic of water than dry 
sand. The drought in the Sahel appeared for a while as a 
figment of tlie imagination. Wot very long after this, it 
was found tliat the data collected refered to the follovjing 
points: 7.3®W, 0°, 7.3‘’E, 17,5°E of longitude and 0° 

of latitude. The three first points very reassuringly were 
found to lie witliin the Gulf of Guinea. A second data set 
was therefore collected and while being more reliable, was 
also more realistic if v/e look at the emissivity distribution 
(Table 5) . Comparing the Figures in Table 5 with an emis- 
sivity of E 0.51, we are then compelled to believe in the 
seriousness of the drought in tlie Sahel region. 


16c .{!) loc (2) lo 

December E = 0.96 E = 0.93 

January E = 0.93 - E 

February E=0.93 E=0.91 E 

April E = 0.94 

June - E == 0.91 E 

Table 5 . 

Computed emissivity for data set II 
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(3) loc (4) loc (5) 

0.94 E = 0.92 

0.95 E = 0.97 

E ^ 0.95 

0.92 


using = E • T^ , 


Conclusion’ 


The model through v?hich we tried to relate brightness 

t 

temperature to moisture content is a N layer stratified 
media with vertical variations in dielectric constant. Soil 
mechanics considerations predict the type of moisture profile 
v?hich governs' these variations. Results are obtained using 
a computer program. 

Even though we did not include in the model any horizon- 
tal variations in the flow equations or temperature variations 
and restricted observations to nadir, the parameters that were 
incorporated in it seemed to be the most important ones as to 
the understanding of remote sensing of moisture in sandlike 
soils. It would then appear that along with being limited to 
a probing depth of twenty centimeters, a limiting factor in 
remote sensing is the sensitivity of the radiometer. The 
results show that unfortunately no appreciable moisture con- 
tent can be determined (while land-sea boundaries have been 
shown to be quite distinct) . However, the mixture formulas 
should be reviewed as well as the correlation betv?een more 
accura-he ground truth data and brightness temperature before 
any final conclusions can be drawn, . 


Data Set I 


The five locations considered are: 


Latitude : 


QO QO QO 


0 ® 0 ® 


Longitude : 


7.3®W 0® 7,3®E 


17.5°E 21°E 


Measurements were all taken between 10:00 am and 

12:30 pm. 
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22 .2 y 

location(2)- ® 
location(5)- x 




V 


290 

280 

270 ■ 
260 . 

250 , 

240 , 

230 

220 

210 

200 •■ 
190 . 

180 

170 ■ 


K 


Q 


(Sf • 






160 i 


150 


20 


1 ..... 1 

1 _ _ . . , ,• . 1 

[ L 

1 ‘ 

'' • 5in ' 

^ An ^ 

n! 


22.2 


31.4 


55*6 

54*9 


location(l )- 0 

looetio„(2). ■ (22Peb-26Pet) 


^ 

60 GHZ 

58,8 
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I .. . • 54.9 

I locaiion(3)- & 

I location(4)- V , ' (22Peb-26l<'eb) 1973 

t location(5)“ 

i 

ii 

ii . • ' 

I 
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location(l)“ @ 
location(2)~ x 


40 


50 


22.2 


50 ^ 

51.4 


60 GHZ 
53.6 58.8 

54.9 


location(5)- ^ 

location(4)“ X ' ( l6Apr-25Apr) 1975 

loeation(5)-’ 





Data Set II 


The five locations considered are: 



(1) 

(2) 

(3) 

(4) 

(5) 

Latitude : 

17 ± 0.5 

17+0.5 

17 ± 0.5 

17+0.5 

17 ± 0.5 

Longitude : 

12 ± 2V7 

4 ± 2W ■ 

4 ± 2E 

12 + 2E 

20 ± 2E 


Measurements were all taken betv/een 10:00 am and 12:30 pm, 
CHANNEL 1 and CHANNEL 2 correspond to 22.2 GHz and 31.4 GHz 
respectively. Even though only months are listed, they were 
all taken in the same period as those in Data Set II . Some- 
times no data was found for our locations of interest. 
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Location (1) 



CHANNEL (1) 

aiANNEL (2) 

December 

281 

281 

January 

279 

279 

February 

280 

279 

April 

290 

■290 

June 

• 295 

295 


Location (2) 

CHANNEL (1) CHANNEL (2) 

Deceniber 277 276 

January - ’ ~ 

February 278 278 

April - ~ 


June 


293 


292 


Location (3) 


December 

CHANNEL (1) 

channe: 

January 

277 

278 

February 

281 

281 

April 

- 

- 

Oune 

298 

298 


Location (4) 



CHANNEL (1) 

aiANNEL 

December 

280 

280 

January 

272 

271 

February 

2 84 

284 

April 

286 

286 

June 
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Location ( 5 ) 


CHANNEL (1) 


CHANNEL (2) 


December - 
January - 
February ' - 

April 285 286 


June 


Frequency 

(GHz) 















' SJOE 

1 

2 

3 

4 

5 5 

6 ' 

7 

8 
9 

10 

11 

12 

13 

14 

15 

16 

17 • 

■18 50 

19 

20 
21 
22 

23 

24 

25 

26 55 

27 

28 

29 

30 300 

31 

32 

33 

•34 777 

35 

36 

37 

38 

39 
.40 

41 

42 

43 525 

44 

45 

46 

47 89 

48 

49 

50 • 

51 

52 

53 ■ 


EOUCiiR. 

KETIL a,HS,CF,FF (100) ,7(100) ,2I,RF{7) 

REikL F (7), EM (7) , FREQ (7} 

COMPLEX ZJ,n,KS (7) ,KH (7) ,RTEM (7) 

COMPLEX K (100) ,EX (100) ,n (100) , 3(99) , Q (100) , CON (100) , CSQRT, 

1 CEXPfCONJG 

format (2F4. 1, 3X,2F5. 3,3X,F4. 1,12) 

READ (5, 5) (KH (JI) , KS (JI) ,FREQ (JI) ,F(JI) ,JI = 1,7) 

DO 111 JI=1,7 

HS=0. 84 
CF=0.26 
2=0.07 
H = 10 

ZJ=(0,0,1.C) 

21=1.0 
S (1 )=Z+CF 
M11=N-1 

DO 50 M=2-,N11 

S (a) =S (1-1) +Z 
CONTINUE 

U= (KU (JI) +? (JI) ) / (KS (JI) +P (JI) ) 

Y (1) = ((H5/5(1)) 

FF(1) = ( (1-EXP (7 (1)) ) -0,1 9) ^0.56 
N11=N-1 

DO 55 I1L=2,Nn 
7 (ML) = ( (H 3/S (ML) ) ^^^ 2 / 2 ) ^ (-1) 

FF (KL) = ( (1-EXP (7 (ML) ) ) -0.19) ^ 0 . 56 

CONTINUE 

DO 300 J=2,N 

JK=J-1 

Q (J)= { (K3 (JI) *U) + (FF (M-JN) * (KN (JI) - (KS (JI) *U) ) ) ) /(U + FF (N-JH) 
1*(1“U)) 

CON TIN U 3 , . ■ . 

K (1) =0. 209*?REQ (JI) 

DO 777 XL=2,N 

K (IL) =0. 209*FREQ(JI) '-i^CSQRT (Q (IL) ) 

CONTINUE 

L=N-1 

S(1)=K(1)-K(2) 

R (1)=R(1)/(K(1) +K (2) ) 

EX (1) = (0.0,0. 0) 

DO 525 1=2, L 
R (I)=K (I) -K (1 + 1) 

E {I)=R(I)/{K(I) +K (1+1) ) 

EX (I) =2. 0*K (I) * (3 (I) -S (1-1) ) 

CONTINUE 


CON (L) = (3.0,0. 0) 

IF (AIMAG (EX (L) ) .GT. 30,0) GO TO 89 
CON (L) =CEXP { (0. 0, 1. 0) *EX (L) ) =i=R (L) 


CONTINUE 

LL=L-1 

IF (LL. EQ. 1) GO TO 378 

JKL=LL-1 

DO 526 JK=1,JKL 

I=L"JK 

CON (I) = (3 ,0,0. 0) 
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54 


IP (AIMAG (BX (I) ) .GI. 30. 0) GO TO 88 

55 


CON Cl) =CSXP (2J=5=EX {1} ) V (H (I) tCOW {I^ i) ) / (ZI+R (I) *CON (1 + 1) ) 

56 

83 

CONTINUE 

57 

526 

CONTINUE 

58 

378 

CON (1) = (R (1) fCON (2) } /(-ZH-R (1) *CON (2) ) 

59 


RTE'1 (JI) =CON(1) 

60 

i 

RF (JI) = (CABS (RTEH (JI) ) ) =^^*2 

61 


(JI) =1-RF (JI) 

62 

111 

CONTINUE 

63 

10 

FORMAT (* *,7F6,2) 

64 


iJRITB(6,10) (SH (JJ) , JJ=1, 7) 

65 


STOP 

66 


END 
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APPENDIX IV 


The Brightness Temperature of a Half-Space Random 
Medium with Nonuniform Temperature Profile* 


I>, Tseng and J. A. Kong 
Department of Electrical Engineering and 
Computer Science and Research Laboratory 
of Electronics 

Massachusetts Institute of Technology 
Cambridge# Massachusetts 02139 


Abstract 


The problem of microwave themal emission from a half-space 
random medium with a nonuniform temperature distribution is solved 
with the radiative transfer approach. For constant absorption 
and scattering coefficients, the brightness temperature is deter- 
mined by a simple close form formula. Physical interpretations 
and numerical results are illustrated and discussed for the various 
cases, \ 


k 

* This work was supported by the NASA Contract 953524 and the 
Joint Services Electronics Program (Contract DAA30 7-74-C-0630) , 


182 


1 , Introduction 


In passive remote sensing of the earth with microwaves, 
the brightness temperature reading of a radiometer depends on 
its angle of observation and on the microwave emissive properties 
of the observed area. Using the model of a half-space dielectric 
n^dium, various theories have been developed for interpretation 
of data collected from satellite and spacecraft. It is well 
appreciated that in snow, ice, or desert areas the subsurface 
temperature profile, absorption, and scattering are dominant fac- 
tors in the surface brightness temperature [Kunzi et al 1975] , 
^Sliming uniform temperature distribution, Gurvich et al [1973] 
derived expressions for the brightness temperatures of ^ half- 
space random medium with a laminar structure. Using a radiative 
transfer approach, England [1974] considered nonuniform tempera- 
ture profile and assumed isotropic scattering phase function. 
Stogryn examined the brightness temperature of a vertically 
structured medium with no scattering [1970] . He also studied 
scattering by random dielectric constant fluctuations in the low 
frequency limit using a perturbational approach [1974] . In this 
paper we propose to solve the problem of microwave thermal emis- 
sion from a half-space medium with a laminar structure and non- 
uniform temperature profile , 


B®^s^0DTTcrBn.m^ of ttu 

fVmaWAI. PACK - 
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II , Formulation 

M 

Consider a half -space random mediiam characterized by a random 
permittivity 





(1) 




(2) 

where e * is a constant and 
m 

e- and e " 
r m 

are 

functions of z 

only. Assume that e„" << 

mm 

and that 

(z) 

denotes the randomly 


fluctuating part whose ensemble average is zero . The medium also 
possesses a nonuniform temperature distribution T(z) . A radiometer 
is sensing at an angle 0^ from nadir. The result of the brightness 
temperature is dependent on tlie polarization and the observation 
angle 8^. Let I^ and I^ denote intensities for the upward and 
downward radiation inside the medium which make angle 6 with the 
nadir (Figure 1) . The angles G and 0^ are related by Snell's 

law. In what follows all results are expressed in terms of 0. 

« 

The radiative transfer equations are written as 


cos 0 ^ = “tc I+tcB_ + »cJ (3) 

e u a T s u. 


d I, 

cos 0 _ = Ja (4) 
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with 


j n -I. I + I 
u 2 • 2 


I , + I„ 
2 ^ 2 


sxibject to the boundary condition 


^ot 


z « 0 


(5) 


z = 0 


where is the Fresnel reflectivity at the boundary which depends 

on 0, 6^, polarization, and the properties of the medium. In (3) 
and (4) 


B « ^ — T(z) 


( 6 ) 


V 7 here T(z) is the temperature distribution, k is the Boltzmann 
constant, and X is the free space wavelength. Equation (3) de- 

i 

scribes the change in radiation intensity as it propagates for 

a distance of dz. It decays by tc^ I^ dz and is reinforced by 
the thermal source ic ‘ B_ dz and the scattering source tc„ J dz . 
A similar interpretation applies to (4) . In (3) and (41 , k is 

C* 

the loss per unit length of intensity due to absorption 


K = 2 k " 
a m 


( 7 ) 
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is the total extinction 


where = Im(k^) « Im /o3^y e^, and 

per mit length caused by absorption and scattering. 


K ~ i: + 

e a 


( 8 ) 


The scattering loss per unit length k and the scattering phase 
fimctions and Pj^ are derived in Appendix A for the TE and 

TM waves; TE waves are polarized with electric field vector 
perpendicular to the plane of incidence and TM waves have polariza 
tion parallel to the plane of incidence. The results are 

(1) • For TE Waves 

k ( 1 + 2k ^Jl^cos^G > 

K = J2 H (9) 

cos 9 I 1 + 4k ^S,^cos^0 ) 

. . ni ' 

1 + 4k ^Jl^cos^G 

P = : - (10) 

1 + 2k ^Jl^cos^G 
m 



1 

1 + 2k ^Jl^cos^e 
la 


** (11 > 


(2) For TM Waves 


REPRODUCIBILITY OF THE 
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k ® Ai (1 + 4k ®j!.^cos^9 + cos^29 
_ n* 1 m 

^ ,2cos 0 1 It 4k„^A^cos^6 


(12) 


C13) 


2 (1 + 4k ^£,^003^6) 

p k; ^ - 

^ 1 + 4k„^£^cos^6 + cos^20 

m 



2cos^20 

■ I ^MMMI - - III I ^ 

1 + 4k„^2^cos^0 ■!• cos^20 
m 


(14) 


In (9) - (14) , i denotes the correlation length of an eKponential 
correlt...ion fmction and & the variance of the flucuation in 
dielectric constant. 

The task now is to solve Equations (3) - (5) for given medium 
properties. Once the upward intensity is obtained, the brightness 
temperature T_ as measured by a radiometer is given by 




T = 


2 k e„‘ 
m 


" ^otJ 


(15) 


z = 0 


For the general case when and are functions of space 

coordinates. Equations (3) - (5) can be cast in the form of inte- 
gral equations and solved by an iterative approach. The formulation 
is presented in Appendix B where a solution to first order in k 
is obtained to check with that derived from exact solutions as pre- 
sented in the following section . 


Scattering Coefficients 


When and of the half-space medium are constants 

a 5 

independent of space coordinates , (3) and (4) become two first 
order differential equations with constant coefficients. We as- 
sume a nonuniform temperature profile of the type 


T(z) “ 


The solution is found to be 


a + K^y cos 0 

I , - P ® •{- B + S 


ct^ - Y^cos^S 


B, e 


“ " *^a ctz/cos e “ ® vz 

I = P ^ + B t B, e^^ (18) 

a + )c^ ° - Y^cos^S 


%} 1 


-ihHH 


P = - 


(a + K ) 

a 

(o + Kg) - rg^(a - Kg) 


- ^ot'> ®o 


(a^ + K V 0) - r . (a^ - k^y cos 9) 
a ' ot a 


■- Y^cos^ 0 
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= ‘ 22 ) 


and the scattering albedo 
2> ~ ‘^s/‘^e’ 


The first terms in (17) - (18) are the homogeneous solution to the 
differential equations and the second and third teirms are the parti 
ctdar solutions . 

The brightness temperature is determined from (15) , 


T = 


2>=a - ^^ot> 


(o + Kg) - r^^(a - Kg) 


T + 
o 


a 


T. 


(24) 


a + ycos 6 


Physical interpretations of this simple result are explored in the 
next section, 

i 
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IV. Numerical Results and Discussions 


We now examine several special cases . (a) For medium with 

uniform temperature distribution and no scattering, Tj^ — 0, 

K == 0 , and a “ ic . Thus 

S M 




(25) 


The emissivity is seen to be given by 1 - This result also 

agrees with that obtained by reciprocity arguments . 

(b) For medium with no scattering, fc„ = u - 0, the result 

s 

is 





'^h 


K + YCOS e 

a ' ^ 


(26) 


As a numerical example, we consider the subsurface ten^erature 
distribution of the Amundsen -Scott station in Antarctic [Lettau, 
19713. The temperature preafile is first fitted with exponentials. 
We find that for December 31st, \ 

T^(z) = 222 + 34 
for August 31st, 

T2(z) 222 - 10 
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and for April 1st, 


✓ 

T3(z) = 222 + 81 - 88 


The unit for temperature is in degree Kelvin. The temperature pro- 
files are shown in Figure 2 , brightness temperatures as a 

function of frequency at nadir are shown in Figure 3 for 

— 1.8(1 + i 0,003)e^. It is seen that at very low frequencies 
Tg is determined by the temperature at greater depth which is 
essentially a constant T^ = 222 with emissivity of ice surface 
approximately equal to e = 0,978. As frequency increases, the 
subsurface temperature becomes more important. At very high fre- 
quencies the brightness temperature approaches the value e(T^ + ^ 

as seen from (26) , 

.It is interesting to note that (24) also reduces to (26) when 
= 0 which corresponds to very high frequency. For when = 2, 

o - K^. 

a 

(c) For medium with scattering and uniform temperature dis- 
tribution, Tj^ “ 0, and we find ** 


T = 


- ^ot> 


"a’ ■ ‘'ot'® - '=a> 


(27) 


Numerical results are presented in Figures 4 and 5 and compared with 
that obtained by Gurvich et al (1973] , We see that for the Shelf 


OF ths 
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glacier [Figure 4] , higher brightness temperatures are predicted 

r 

with our model. This is because their results arc applicable to 
low scattering cases while the scattering albedo for the Shelf 
galcier is rather high. In the case of continental glacier [Figure 
5] where the scattering albedo is quite small, both models agree 
^uite well except at resonance where the scattering is largest 

(d) In the case of uniform temperature distribution and small 
scattering albedo, we expand (27) to first order in u by noting 
that 

a ^ Kg(l - to Pf/2) (28) 

^d 

“ K^(l *- w) . (29) 

The result is casted in a form to compare with that obtained by 
Gurvich et al [19731 . 


8k^‘* cos 6 (1 + 4k ^fi,^cos^0) 
m m 



This expression differed from their result in which there is a 

cos 0^ in the numerator instead of a cos 6 in the denominator 
o 

as shown in (30) . An alternative derivation by using the v/ave 
approach is given in Appendix C which confirms (30) . 
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*' / ' 

Numerical results are given in Figure 6 to compare TE and 
waves as a function of radiometer viewing angle 0^. Clearly# 
scattering lowers the brightness temperature at all viewing angles . 

For non uni form temperature distribution# expanding (24) to 
first order in w yields- 


Tg = (1 


^ot> 




1 - 




T. 


h 


1 + ycos 0 /k 


1 “ 


ycos 


TS 

w P, 


1 + ycos S/k 


w P, 


4 


(1 - 




(31) 


This result is checked with the iteration approach as presented in 
2ippendix B , 


(e) Finally, we plot Tg for the temperature profile for 
the Amundsen -Scott Station by using the exact formula in (24), The 
results are shown in Figure 7, The radiometer is viewing from nadir. 
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Conclusions 


The problem of a half-^space random medium with nonuniform 
temperature distribution has been solved with the radiative trans- 
fer approach. For constant absorption and scattering coefficients, 
the result in (24) is very simple and admit direct physical inter- 
pretation. It is interesting to note by comparing (24) and (26) 
that a corresponds to an effective absorption coefficient in the 
presence of scattering. Various special cases are discussed. For 
medium with inhomogeneous absorption and scattering coefficients, 
an iterative approach can be used. In Appendix B, such an approach 
is illustrated and shown to give identical results in the limit of 
constant k and k . The theory and results are compared with 
previous approaches. To check the results that we obtained, the 
wave approach is used (Appendix C) . It is a straightforward matter 
to extend the radiative transfer approach to multi-layer models. 


* 
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appendix A 


In t±iis appendix, we derive the scattering loss per imit 
length k and the forward and backward scattering phase functions 
and as given in (9) - (14) . The scattered electric field 

E® of an incident field E^ can'be expressed in ternis of total 
field E by the wave equation 


V >c V X E Cr) - kj^^ E®(r) « e^(r) E(r) . 


(a-1) 


In terms of the dyadic Green 's function we have 


E®(r) = G{r, r') • e^frM E(r*) 


(A-2) 


Inside the integrand in (A-2) we write 


E(r') E^ e 


ik • r 


(A-3) 


where k is the incident wave vector, and A is a unit vector 
denoting the polarization of the incident wave. Using far-field 
approximation, we obtain 


r<5) = (1 - k 1 

s s o 


4irr J 


d^r* e„(r'). e 


i(k - k„) • r’ 


(A-4) 


and 
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m 


o 






sin^X 


16 rr^e„^r^ J 
m •' 




i (k - k_) • (r, - r,) 
e ® i. ^ ^ ^^_ 5 j 


where x angle between the incident polarization and 

A ^ mm « 

the scattered direction k^. Letting p - r^^ - r 2 and 


^ ^^^1 ” ^ 2 ^ (A- 6 ) 

we find the bistatic scattering cross section per unit volume to 
be 

. k ^ A 

q(k» k ) = sin^X — r 

® 16 TT^ 

This, result' checks with that obtained by Gurvich (1973) and Tatarski 
(1961) when we identify A “ 4 ct^q^ • ' 

To calculate scattering loss per unit length, we consider the 
case where 

b(p) = . (A-8) 

i 

A X 

.Let the plane of incidence be the x - z plane. We write 

A A 

k = X k sin 9 + z k„ cos 0 (A-9) 

m m 

k = X k sin 0 cos c()_ + y k sin 6 sin + z k cos 0 . 


/ 


d®p b(p) e 


i(k - 


J^s> 


(A-7) 



The scattering loss per unit length is 


/' 




, kg) 


k^^A 


2Z 


/ 2tT rTT 

d4>s / sin e sin^x — 

Q / Q it 1 2 ^ i*Ti^f? CS MB A \ 2 (J ^ 


4 1 + k„ (cos 6 “ cos 

m s 


6 (kj^ sin ^ sill Qg cos (J)g) 6 (k^^ sin 6^ sin , 


(A-11) 


Within the range of integration, the Dirac delta functions give non- 
zero results at (1) ‘ ^ =0 and 6=6 which corresponds to forward 

scattering and (2) = 0 and 0 = tr - 0 which corresponds to 

S W 

backward scattering. This finding is essentially dictated by phase 
matching for the laminar structure [Kong, 1975] , 

The integral in (A-11) is evaluated by changing variables 


u = k„ sin 0^ cos 4 - k„ sin 6 
m s ^s in. 


(A-12) 


V = k sin 0„ sin 4 
m s ’^s 


(A-13) 


which gives the Jacobian l/[kj^^sin 0^ cos e^j , 
(1) For TE waves, x “ Tr/2 , and 


hj LI 

K_ = — —a + 

2 cos 0 


1 + 4k„^£.^cos^0 
m 


1 


(A-14) 
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The first term is from forward scattering and the second term from 

'backward scattering. Equation { 9 ) follows from (A-14) , Equations 

(10) - (11) are obtained from (A-14) by noting that = 2 

and P^/Pu = 1 + 4 k„^)l^cos^0. 
f b m 

(2) For TM waves, we have for forward scattering, x - it/ 2 and 
for backward scattering sin^X " cos^20. The scattering loss per 
unit length is calculated to be 


K 


s 


2 cos 6 


cos® 2 9 

1 + - . , 

1 + 4k„®^®cos®0 
m 


(A-15) 


Again the first' term, is from forward scattering and the second term 
from backward scattering. Equations (12) - (14) are direct con- 
sequences of (A-15) . 
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Appendix B - 

t 

In general, Equations (3) to (5) can be solved by an iteration 
approach. First we cast the Equations into integral ■ form 


dz ^ (C ' fz ' ) 

e /• z 


I„(z) = a 


+ I ^ dz ' (k^(z‘)B^{z*) + tc^(z') J^(z’)) 


■/ 


z* 


dz” K’ (z“) 
0 


(B-1) 


■/: 


dz“ X* (z") 

iLa 


I^(z) -J dz ' (k^(z')B^(z ') + Kg{z')Jj^(z') e 


r 


+ ^ot / ' (tc^<z ')B^(z*) + K^(z')J^(z') 

i«CO 


-J* dz" fc^(z") -y" dz" 


K' (Z'V) 
0 


(B-2) 


where k„ ' = ic_/cos 9, tc^' = k_/cos 0, ' = x„/cos 0. For 

O 0 cl 3 0 0 

iteration we assume ' is small and write 

s 


Iu(2) = + ... 


(B-3) 


+ li^' (Z) + iP’(z) + ... . 


(B-4) 
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Substituting in (B-1) to (B-2) and equating the appropriate orders. 


MQ obtain, for m=0, 1, 2, ••• 


(m + 1 ) 


dz' ' izn e 


dz“ K*(z*') 
£ z* ® 


(B-5) 


dz’ K^(z’) (z 


/ 
•) e = 


dz" tc'(2") 


+ ^otf^ ( 2 ') 


dz" k'Cz 


") -p 

J 5 - I 


(B"6) 


®o ® 


-/. 


dz* K‘{z') 
>»• 


dz ' K ' (z') B„(z ') e 

* I fC 


dz" k’ ( 2 ") 


(B-7) 


dz * K^(z ')By (z ') e 


dz" K* (z") 
6 


t r^^ \ dz' k^(z']B^(z') e 


-/ 


dz" Kl(z") 
■ 0 


dz" K^(z") 


(B-8) 


Consider the special case v/hen k' and jc ' are independent 


of z. Assume a nonuniform temperature profile of the type 


T(2) = e^®. 


CB-3) 


We find from (B-7) 


(z) 


kJ k i + y 

e e ' 


(B-10) 


and are given by (21) and (22) . 


From (B-8) , 




yz , 

+ e ' + 

K ' - Y 

e e ' 


'k' kI b. 
a o . an 

■■- ■ ■ ' -t* -■ — ■ 


Kg “ Y. 


. _ , =o , *'1 ®h 

+ [ — r~ 


<1 Y. 


K'Z 

e 


(B-11) 


Introducing (B-10) and (B-11) in (B-5) , we obtain 


( 1 ) ! 


u 


W tc * B_ K ' 

a o ^ s f 


''a ®h 


z = 0 


+ Y)^ 


S P _ 
2 K 


a n 

• 2 „ „3 


tu P, 


' B K ' P 
_a o ^ a n 


■' 


+ r, 


ot 


^ k 1 B^ fC* B, 
a ^ a h 




. (B-12) 


The brightness temperature is then given by 




XMl - 


^B = 




2ice^‘ 

m 


( 1 ^“’ ( 0 ) 


+ I 


( 1 ) 


U 


v?hich gives (31) . 


( 0 )) 


i 
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Appendix C 


For uniform temperature distribution with constant the 

brightness temperature is given by 


■ 'ot ■ 

where s denotes contributions from scattering. Consider a 
plane TE wave incident upon the half-space medium at the in- 
cident angle 0 , 

i (k X ^ k z) 

' <C-2> 

The equations that govern the electric fields in the upper and 
lower regions are 



0 




= Ef(z) 


« 


i 


(C-3) 


(C-4) 


where k^^ + = k^ / ~ ^m^ ' t electric 

field in the half -space medium. 

We first solve the Green's function for (C-3) and (C-4) and 
then use Boam approximation -to calculate the scattered fields. 
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The Green's functions in the two regions are goveimed by 



(C-5) 


(C-6) 


Imposing the following boimdary conditions at z « 0 : 


ot tt 


(C-7) 


* ®ot _ °tt 


dz 


dz 


(C-8) 


and the radiation condition at z we find 


1 - R . ik z - ik z ' 
G^^(z, 3') = 2t.e 2 


2i k. 


(C-9) 


mz 


where is the Fresnel TE reflection coe^fficient between two media 

with permittivities c and e^. The Green's function in region t 
G^^(z, z*) is given by a more complicated expression. 

The Bom approximation to the scattered field in region 0 is 





dz 


G (z. 


z ') e-(z') E (z*) . 


(C-IO) 



The zeroth order field amplitude in region t is given by 


-ik ' 

V'*’’ = ^ 


(c-ii) 


Using the correlation function in (A-a) r we obtain 


^ (1 - Ak 

<[e®I2> E 2 


mz 


/>-j: 


dz. 


0 0 


E, 




1 + 

mz 


r\j E 

O 


(1 - r^,)^ A k ^ A 
2 ot m 


8k^"cos 0(1 + 


CC-12) 


in view of the fact that ■ k^^" = k^' k^.% k^- v k„ and 

k ' '>> k 'i' k cos 0 . Equation {C-12 ) gives the second term in 
mz ^ mz ^ m 

(30) . 
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Figure 1. 


Figure 2. 


Figure 3. 


Figure 4. 


Figure 5 . 


Figure 6 , 


Figure 7. 


Figure Captions 
Geometry of the problem. 

Subsurface temperature distributions at the Amundsen- 
Scott Station in Antarctic. 

Brightness temperature (without' scattering) as a function 
of frequency for the three temperature distributions 
shown in Figure 2 , 

Brightness temperature for shelf glacier and compared 
with that obtained by Gurvich et al and depicted in 
broken lines . 

Brightness temperature for Continental Glacier and 
compared with that obtained by Gurvich et al and de- 
picted in broken lines . 

Brightness temperature as a function of viewing angle 

, i 

for TE and TM v/aves . 

Brightness temperature as a function of frequency for 
different scattering media for the three different 
temperature distributions in Figure 2 . 
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Microa^ave Remote Sensing of a Two-Layer Random Medium^ 


hi 
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Computer Science and Research Laboratory 
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Abstract 

Modified radiative transfer (MRT) equations are derived 
for a tv;o-layer random medium from the Dyson equation by using 
the nonlinear approximation and from the Bethe-Salpeter equation 
by using the ladder approximation. The MRT equations yield 
simple and useful solutions which are applicable to both active 
and passive microwave remote sensing. 



I * In troduction 


In microwave remote sensing of the earth or other planets, 
theoretical models are essential to the interpretation of the col- 
lected data. ®ie effects of absorption and scattering have long 
been recognized as dominant factors in both active . and passive 
sensing. Using the model of a half-space random medium with laminar 
structure, Gurvich, Kalinin and Matveyev^^^ investigated thermal 
radio emission in Antarctic areas. The problem of scattering by a 

medium with a small random fluctuating part in permittivity was 

f21 ■ 

studied by Stogryn . Assiuning a half-space containing random 

[31 

distributed scatterers embedded in a low loss dielectric, England 
examined emission darkening caused by scatterers with a radiative 
transfer approach. 

The problem of random media can be treated with the radiative 
transfer method which deals with energy fluxes and with the re- 
normalization method which deals directly with field quantities , 

The renormalization method gives rise to the Dyson equation for the 

mean field and the Bethe-Salpeter equation for the covariance of 
[3 41 

the field, ' ■' In solving the equations, the bilocal approximation 
is usually applied to the Dyson equation which is then solved by 
mathematical techniques such as the Fourier transform method. A 
ladder approximation is made on the Bethe-Salpeter equation which 
is then solved by the method of iteration. In the case of 
multiple wave scattering, the method of iteration involves solving 
many integrals and leads to complicated results after one of two 


iterations . Under the assumption of far field interaction and 
incoherence among waves in different directions/ radiative trans- 
fer equations have been derived from Be the-Sal peter equations to 
study multiple scattering. 

In this paper we investigate the problem of scattering and 

emission of microwaves by a slab random medium with a laminar 

structure and bounded by different dielectric on each side. We 

employ the nonlinear approximation rather than the more 

popular bi local approximation to the Dyson equation. A two vari- 

ri51 

able expansion technique is applied to obtain the mean Green's 
function, which is then used to derive a set of modified radiative 
transfer (MRT) equations. They are modified because in our deriva- 
tion we include coherent effects which are important in the case 
of bounded medium. The MRT equations are then solved and illustrated. 


II , Formulation 


Consider a tv/o-layer medium with boundaries at 2 = 0 and 
2 = -d (Figure 1) , The slab between the two boundaries is a random 
medium x-^ith a permittivity 

(X) 

<c^^iz)> = 0 . ( 2 ) 

The mean permittivity is independent of depth. The random part 

of the permittivity is real and has a sero ensemble average 

as indicated by (2) . The variance of the fluctuating part of the 
permittivity is also assumed to be small, 

A monochromatic plane wave is normally incident upon the two 
layer medium 

= exp(-i k^z) , (3) 

The time dependent factor expt-iwt) has been suppressed. The 
governing equations are the one-dimensional wave equations in each 
region . 



in region 0 , 
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<EE*> 


2 


+ < 86 *> 


field intensity is 


The Dyson equation can be obtained for the mean field in region 

n[4-53 


d^ 

da^ 




where Q(z/ is the mass operator v/hich is an infinite series 

of the unperturbed Green's functions Zj) and the correla- 
tion functions a^) = <f{z) f(z^)>. Writing in terns of the 

mean Green's function 2^), the Dyson equation becomes 




'11m 


r ® 

(Zf z’) = / 

•^-d 


da, Q{Zf 


Zi) 


^llm^^l' 


a') 


+ 6 ( a - a ' ) . 


( 11 ) 


In writing the Green's function, we use the first subscript to de- 
note the region in which the observation point is located and the 
second subscript to denote the region in which the source is located 
For the covariance, the Bethe-Salpeter equations takes the form 


<f^(z) = J / dZ2 / ' / dZj- Zj^) 


^llm^^'" ^ 1 '^ ^ 2 ' ' 


^2 ^^lm^^2^ 


,* 

'Im' 


E. ,(z') 




(12) 
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. Thus 


where ^{2^, z^r 2^', is the intensity operator 

we see that the Dyson equation must be solved before one can solve 
the Bethe~Sal peter equation. 

Both the Dyson equation and the Bethe-Salpeter equation are 
exact equations. However, the mass operator and the intensity 
operator are in the form of infinite series. Two types of approxi- 
mations, the bilocal approximation and the nonlinear approxi- 

[13-14] 

matron are used in solving these equations . In the folloviring 

sections, we apply the nonlinear approximation. As we shall see, 
this approximation preserves the principle of conservation of energy 
and enables us to derive the modified radiative transfer equations 
for bounded medium. 
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Ill . Mean Green’s Function 


In order to solve for the mean Green's function, we impose 
the nonlinear approximation on the mass operator in (11) 


Q(z, Z|j^) Qj^ ( z , — ^llm^^^ ^1^ • 


(13) 


The correlation function ®^(z, z^-) is assumed to take the form 


Bj(z, z^) = <fCzl £(Zj^)> = « e 


-! z - 


(14) 


where 6 is the variance and £ is the correlation length of the 
fluctuating part of the permittivity. 

Before proceeding to the solutions, we define soite parameters 
and specify their ranges of interest to us. Let k^^^ ^ + i k^^^ 

with k^^ and k£^^ both real (we use prime and double prime to 
denote the real and imaginary parts of a quantity respectively) , we 
require k^^ << The absorption v;ave nun±>er and the absorp- 
tion length £ are defined to be . 

a. 



(15a) 


(15b) 


The scattering wave nuv.iber k and the scattering length £ are 

s s 

defined to be 
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(16a) 


*c_ = 




= 


2TT/Kg . 


The v;avelength inside medium 1 is taken to be ~ 2Tr/k^j^^, 

The range of interest to us is for << 

The Dyson equation under the nonlinear approximation becomes 






-| 2 - Zil A 

® - 2‘) . 


(17) 


This equation can be solved by using the two variable expansion 
f 15 1 

procedures. •' V?e define long distance scales K ~ fiZ/ = 6z’, 
and ^^om ^ = 62 it follcscs that 


_d^ 

dz^ 


-Jl + 26 + 62 

3z2 3235 



(18) 


The mean Green's function is expanded in perturbation series 






C‘) = 


'^llmo^^' 


e; z*. Kn 


+ 6 5’) + ... . (19) 
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Substituting (18) and (19) in (17) , we find that on balancing 


terms to zeroth order in 5, 




and to first order in 6 


T ‘ ^"im '"11ml 


+ 1 (z, e? z', r) = -2 _L_.g^^^^(z, Kj 2 


3z3C 


+ I ^ z,) e ' 1'' Gii„„(Zi, z 


-| z - Z I/A 


^Im ""llmo'"" "-I- 


’llmo^n 


The solution to (20) is v/ritten as 


[U(0 +W(C) e ] (22 


for z > z ' and 


tU(E) + W(E) e 
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for z < z'. In (22), note the syiranetry property of the Green's 
function. 

Following the two variable expansion procedure, we substitute 
(22) in (21) and eliminate secular terms. We observe that , 

B(i^^), U(5^) and W(^^) as variables in the integrand are varying 

on the i and % scale. They vary much more slowly than 
exp(“|z - Zj\/l) V7hich is 1 at z ~ and decays on the I 

scale. These terms can be taJ^en out of the integral sign and siab- 
stituted with their values at We find that in order to 

eliminate the secular terms, the follov;ing equations must be satis- 
fied: ■ 


g MO 


2 + i3k, Jl\ fl - ikn^il 

B(0U(5) j + A(OW(?) ' 


1 t i2k3L^£ 




= 0 


(23) 


aB(g) 


- g B(e) 


2 - i3k, fi,\ /I + 

A(?)w(e) I ±HL)+ b(c)d (0 ' 


1 - "2k,/ 


"2k,/ 


0 


(24) 


g U(0 


1 + ik, £ \ /2 - i3k, £ 

B(e)U(5) { ±5L_]+ A(5)W(?) ' 


^ "2^im^ 


1 - i2k,/ 


= 0 


(25) 


dW(0 


~ g w(5) 


1 - ik, £ \ /2 + i3k, £ 

A(OW(C) [ — — + B-(^)U(?) ' 


1 - i2k,/ 


^ "2k,/ 


= 0 
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(26) 


where 

g = (27) 

Combining (2 3) with (26) , and (24) with (25) , we obtain 

W(^) + A(^) = 0 ' (28) 

de • . 

and • ■ 

B(^)‘ + u(^) - 0. (29) 

d^ . , 

Thus 

A(?)W(C) = 

B(e)U(C) = L 2 

where and L 2 are constants independent of ^ . Substituting 

(30) -(31) in (23) and (24), we can solve for A(^) , B(^) , U(^), 

and W(^) . 

-g{M-+M^)e 

MK) = . (32) 

g(W- + N-)? 

B(5) = C 2 e - (33) 


(30) 

(31) 
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U{?) e ^ 


L- g(M + 

W(5) - ^ e 


where 


“i = 






L, 

^ d 


2 + i3k^^Jl 


1 + 


^1 == ^1 






^2 = 4 


^ ■** 




Substituting in (22) yields 


(34) 


(35) 


(36) 


(37) 


(38) 


(39] 




^ 9(^2 ■*■ M^)<Sz ^2)62 “^^1] 


+ e 


^g{N, + N„)6z' ik,„ 2 ' 


L- g(Mp + M,)6z 



far z > z • and 


z') = {C e 


-g£M2 + M^)6z' g(N^ + N2)6z' “ik^^z*^ 


-g(N +N2 )«z iki„= g(M2 + H^) 6 z 

{L* e e + — e el 

^ C 

(41) 

for z < z', where C = C^/C^ in (40) and (41) is a constant to 

be determined together with and by the boimdary conditions . 

The Green’s functions in (40) and (41) are continuous at z = z'. 

By (20) we must also have z = a’ " ^‘^^llmo^‘^^'^ z = z' “ ^ 

+ — 

which gives 

- L 2 = l/i2k^^. (42) 


We next match the continuity of z = 0 

and at z = -d. At z = 0, we have 

where is the reflection coefficient at z = 0 , At z = -d 

we have 


C 1*2 g[N^ + + M2l6d 


' i2k, d 
= ^2 -■ 


C44) 
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where R 


is the reflection coefficient at z = -d 


12 


p — _ Q — — p * 

k.„ + k ■ 01' 

Im o 


- k, 

V = -IB 1 

■'im + ■'a 


We de fine 


''l ' ^ "l" “ " ’'ira ■'■ ^*2^* 

Ii2 H na' + i Ha" = ki„Ci - 6MNj^ + n^)). 

In view of (36) -(39) and (42)-(44)| V 7 e find 



(4! 

(41 

(4( 

(4( 

(4’ 


ho ^^lni ^ 


‘Im 


3 - i4k, £ 
lin 

1 - i2kj^„a 


4- R, 


10 


^2 


3 + ^4k^j^S, \ i + n^) d 


1 + i2k^^Jt. 


Prom (4 2) and (4 4) , v;e find 


BEPRODUCIBILrry OF THE 
ORIGINAL PAGE IS POOB 


= 


^12 


i(nr+ n-,)d ■ ‘ 


(4‘ 


e 



Thus the zeroth order mean Green 's function as determined from (40) 
and (41) reads 




10 


-in-iz: * i'l-i 

' («12 ® 


i (n, + 


-a 


+ e 


(50) 


for z > z* and 


^llmo 


i^pZ* 

(Zf z') = L^(e 


+ R 


10 


“irin z' i tn.T + 

e ) { 2 ® 


tl2) ■*' 


irii z 

+ e 


(51) 

for z < z'. From (47) and (48) v/e see that the effective propaga- 
tion constants and rt 2 depend on the boundary conditions as 

well as the properties of the media. 

The two variable expansions technique can similarly be applied 

to the Dyson equation for the mean field 


{■£ * ‘i-1 


r° 

+ k?„ I E,„(z) = dzj^ 




-|z - Zj^l/J 


%m<==l> 


(52) 


under the nonlinear approximation (13) . The zeroth order solution 



^Imo 


(z) 




~ir\2^ + ^2^^ *** 

te +^ 12 ® > 


(53) 

where = 1 + Rq^. We also define = 1 + R^q# where 

^10 ^ ^ 01 * 
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IV« Modified Radiative Transfer (MRY) Equations 

We now proceed to derive modified radiative transfer equations 
from the Bethe-Salpeter equation (12) , 

We apply the ladder approximation to the Bethe-Salpeter equa- 
tion (12) . 


l(z^, Zj; Zg') Bj(Zj^, Zj^-) - Zj) - Z2') 


(54) 


Under the ladder approximation (54) , (12) becomes 


enz')> = ski^ /° dzj^ /“ dzg g^3^^(z, z^) g*^^(z', zj) 


-|zi - Z 2 I A 




1X^1/ 


(55) 


The fluctuating part of the field is written as a summation 
of upward and downward propagating waves, 

<$^(z) = e + <£' 3 _^j(z) e , (56) 

We then have 
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ik,' (z - z') 

<f^<z) = J„(Z, Z') e Z') e 


u 


ik' (z+z*) -ik' (Z + Z’} 

+ J ( 2 , z'} e + J (z, z') e 

Cl ^2 


(57) 


where 

J^(z, z') = <fj^„(z) g;L^*(z')> (B 8 a) 

■J^(z, z') = <£j^a( 2 ) ^^^*(z')> (58b) 


J (z, z') = < 63 ,,j{z) fj^^*(z')> (58c) 


(z, z‘) = <fj^^(z) (58d) 

2 


J and J represent the correlation between the upward and 
downward propagating waves , 

We wish to derive radiative transfer equations that govern the 
behavior of intensities 


1^(2) — J^(2/ Z) 
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I^(z) = J^(z, z) 

• 

(60) 

I (z) = J (z, z) 

= J* (z, z) . 

On 

(61) 
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^ Ke restrict the separation of z and z’ to be on the wave- 
length scale which is much lass than and scale, jz - z'| 

also require that |z|, |z'|/ | 2 + d[ , 


z ' + dl , d >> A 




S,f namely both z and z' are far away 


from the two boundaries on the wavelejigth scale and the correlation 
lengtli scale. However, both z and z* can be close to tlie boun- 
darii- on the scale, and so are J„, J,, J and , 

d 5 d d G ^ 

without loss of generality we consider (55) for z > z*, balance 
terms by their phase factors and then let z' z. After consider-. 
able algebra (Appendix I) , we find 


-K„z “(k^ + K^)d - K , 2 

I^(z) = e ^ l^izr z) + rj ^2 ® 


(62) 


K, Z 


K^Z 


I^(z) " e ^ Ij^(z, z) + e ^ ^ 4 t 2 , z) 


(63) 


i2k,‘ z 

I^(z) e = -R* 3 , e 


Krij + ni*)z 




whe re 


i(r|, +119)^ + i(ni + Ti9*}z 
+ ® 14 ( 2 , z) (64) 


K_ K D (3 -1* 8k^ 5,=^) 
Ki ^ H. -S + s 1 Im 


2 2 1d|=^ (1 + 


(65) 


Di (3 + 8kJ 

« 2ri2'' = ® ^ 


2 2 |d|^ (1 + 


( 66 ) 


■ ° ^ ®01 «12 ® ^ ^ 


(67a) 
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I (z, z) = i 2 / dz, 


I^nll" 


01' V ="^"1 * ''2’^ **=1 ■*■ "2’®! 

— '3b r^2 ® 


D 


-(k. + K«)d (ic, -s* !C^)3, 

+ gfd + xj^2 ® > + 3b >^01 ® 


3f(Koi %2 


iCri-j + n^)d -Krii* + nn*)^ 

e ^ + R* R* e •*■ ) ] 


fC 2 

■'■ ^ 01 '3f * 3b la'^l>> ® ^ 


3 f **01 


lo<==l> 




- 0 


- 32 *) 


- < 3 , 


R* 

^01 




i(n. 


n^* - 




(69) 
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1 + 
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(70) 


(71) 
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and ^12 ” From (64), we can express 

in terms of I„ and I,. 

u a 


i'24„ 


“ r\. 


nj,*)z 


Dn 


£- Hg3_ I^(.) 


k:2X 

e 


" ^2 ^01 ® 


i(Tli ^2^*^ *** ^ 2 ^ 


^01 ^12 ® 


- (k^ + < 2 ) d - K^Z 


i CriT + Ho) “ K, z 

+ Rj^ 2 ^d ® ^ ^ . (72) 


We note that in view of (62)-(63) and (68)-(69), 

■ lafo) “ =^01 iu<o> 
i„(-a) = r^2 ^a*"^' • 

These are the boundary conditions satified by the intensities 
and r^. 

Writing the integral equations (62) -(63) in differential form, 
we obtain 
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dl. 


u 


K. 


dz 




|D| = 


*s Pb '=2^ ^ 

e + 


=sPf 


Pf + ) ^12 ® 


“(iCj^ 4- iC2)d - K^2 


(74) 


dl , K. 

— =^e -— (Pf^a-^PbV 


|T, 


01' 


D 


Pb “(*^1 + ** 


12 


^sPf\ ^2^ 
* I Pf - -^J e 


(75) 


where 


s 


1 + 


(76) 


K « k 

^ iDl^ ^ 


(77) 




(78) 
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(79) 


Pf 


1 


Pb = 


1 + 


Equations (74) and (75) together with the boundary conditions in 
(73) are the modified radiative transfer equations derived from 
the Bethe-Salpeter equation under the ladder approximation. 

In the case of a half-space random medium*- we let d ^ and 
the modified radiative transfer equations reduce to 


dz ® “ 2 


f Pb ^d> 


Pb 


’’oil' ® 


(76) 


= *^0 ^d “ T tPf ^d Pb ^u^ - 


K z 

’’oil' ® 


(77) 


We remark that the principle of energy conservation is perserved 
with the use of the nonlinear approximation for the mean Green's 
function instead of using the more popular bi local approximation. 
Ihis is shown in Appendix II. 
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V. Solution to the MRT Equations 


The modified radiative transfer equations (74) and (75) are 
simply two simultaneous first order differential equations with 
constant coefficients. They are readily solved and written in the 
following form; . 


|T i|2 -a(z + 2d) 

I « -P-f— (P e“^ + Q f„ e . - e ^ ) 

^ 1d|^ 2 


. (78) 


= — ■ (P e + Q e - r-j ^2 ® 


-Kj^Z - (»Cj^ + K2)<3 


^ |D| 


where 


(79) 


P = 


(1 - r3_2)D^ 


^2 ^ 01 ^ “ ^2 ^ 12 ^ “ ^^2 “ ^ 01 ^ ^^2 “ ^ 12 ^ ® 


-2ad 


(80) 


Q = -P 


^2 “ ^12 
^ “ ^2 ^12 
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(81) 


^2 = 


2QC3 - cc) -K3 

■ ^s Pb 


(82) 


'Vi/2,. w . w .1/2 
a - (1 “ to) ^ (1 + - pj^) ' 

2 2 


(83) 
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and (0 “ /k is the albedo • 

s e . 

To calculate eniissivity e of the tv/o layer medium, we note 
that by reciprocity 

e = 1 - r (84) 

where r is the reflectivity. The field in region 0 can be written 
as 



"ik z 


ik z 



E - e ° 

o 

t r 

e ^ . 

(85) 

The 

reflectivity 

•is 

r = <rr*>. We observe 

that <E E_*> = <E,E *> 
o o 11 

and 

= 

jEi^(0)[^ at the boundary 

z = 0 . By (53) and 

(85) 

, we have 


i(rtT + r\'jd 



R + 

01 

<r> - 

®12 

e ^ ^ 

(86a) 


D 


= 1 + <r> t <r*> + <rr*>. (S6b) 

By (57) and (62) -(64) we have 


<EiEi*> 

z 


= |E, (O)l' + 

= 0 ^ z == 0 

- 1 + <r> + <r-^> + l<r>|^ + i,,(o) + i,(o) + i^co) + i *(o) 

Vx d C O 

“ 1 + <r> + <r*> + |<r>|=^ +*i^(0) |t^q[ ^ , (87) 
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Equating (86) and (87) we obtain 


r = <rr*> l<r>|^ + i^(o) |t^q|^ = 


^01 *** ^12 ® 


i(nn -J- Tlo)^ 


D 


|Di=' 


(1 ^ 2 ^ 12 ^ ^^2 " ^ 


“(jCl + K: 2 )d 


12 


) - (f^ - 


• 2*"12 


“(jc^ + _2cta 

i ) e 


(1 - (1 - f^r^,) - (f^ - rm)(f^ - r) 


'2 12 


01 ' 


-2ad 


(88a) 


where 


^01 ^ “ ^*01’ 


(88b) 


We now consider several special cases; 

1. A half space random medium. We let d and (88) becomes 


r = roi + 


^01^2 
^ ■ =^01^2 


(89) 


which is identical to that derived from a phenomenological radiative 
transfer approach. 
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2, The second layer is perfect conductor. We have = -1 and 

^12 ^ Equation (88) becoittes 


r = 


Rqi “ e 


i{Ti:L ’^2' 


1 - e 




“(k, + KT«)d -2ad 

^01 ^^2 " ® ^ + (1 - e ^ ^ ) e } 

— 2ad 

{1 - fg^Ol ® Hi - 1^01 ® 


^^1 ^ Ti2)d 


(90) 


It is interesting to note that as w 1, we have a *> 0, f 2 1 
and from (90) r- ->■ 1 . Thus all the incident power is reflected. 
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VI, Conclusions 


In this paper modified radiative ti'ansfer equations arfs 
derived from the Bethe-Salpeter equations. The nonlinear 
approximation is applied to the Dyson equation and a two variable 
expansion technique is used to find the mean Green's function. 

The modified radiative transfer equations accomodate coherent 
effects which are essential to bounded media and appear in 
observed data, Tlie MRT equations are solved and the results are 
presented in Equation (88) , In Figure 2 and 3 we illustrate 
the emissivity for a half-space random medium and a slab random 
medium with 6 = 0,1 and 0.2 . We see that in the half-space 
case the null in emissivity due to scattering occurs at 
In the slab case scattering dampens the interference pattern and 
decreases emissivity in general. The existence of the interference 
patterns depends on the location of the si±>surface and the extinction 
loss of the random medium. 


Appendix I 


»liODUCIBILlTY OF THE 
PAGE IS POOR 


On substituting equations {50) and (51) into the right hand 
side of equation (55) , we find that 


jgi , • , irio^ ” iri2*2' “in^z + iTi,*z* 

<e 2 ^(z) 6j^*(z')> - {e ^ ^ j. V. « 1 '1 


+ e 


in^z + in, *2' -iri.z in^^z' 

- RJl e 2 1 _ ^ '1 >12 


> l 3 ^(z, z') 


{Rf2 


+ ri2*)<3 + iri2Z ^ in,*z 


- RqI e 


“-iri^^z + in2*z* iri2Z + in^ 
+ e 


I 


-iCn,* + - irinZ - irti*z' 

Ro 1®12 ® ^ ^ } Ijtz, z') 


+ {K„ 3 , e 


+ ^2^^ "** “ in 2 *z' -iTi^z + in_*z' 

- R* Z 1 

Kqi e 


“^^2^ “ iri 2 *z' iCn-, + ri 2 )fi + in-[Z + in *z» 

^01^12 ® > 


T / . r "^'^1 *^ 2 ^^ ■** “ iriT*z' 

l 3 (z, z') + {r ^2 e ^ ^ ^ 




-iri^z + iri2*z* 

+ e 


+ R, 


12 


i(rij^ + H2)<1 + inj^z + il2*z‘ 


-i(n,* + n,*)'3 - in,® - ^ , ,, 

+ RJ 2 e ^ ^ ) I 4 U, z') 


(A.l) 


where 


X^(Zf z’) “ 


« Hm 


41d 12 •'-d 


f dz M dZrt [Rq Q s 


i(tl^ + r]^)^ + ini^i ^ 






+ e ] exp - 


“ ^2 




(A. 2) 


I^Cz, zM - 


iH f dz, r dz, {e 2 ^ - Rg, e ^ 

41dP -^-d z' ^ 


(Ri 2 d 


i(n, + r 2 >'^ _ L 1^1 ' " 2 I 


+ e ) exp ~ 




{A. 3) 
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T^(Zf 2*) 


6 k? .0 z* iri„z, 

iHf dz, f dz,(e 2 ^ 

4|dP ^ z ^ ■'-d ^ 


®01 ® 




{R£2 e 




) exp - 




(A. 4) 


5 k? 0 rO ixio^T “iriTZi 

I.{z, z‘) = dz^J dz (e ^ - R e -^ -^] 

^ 4|D|" ^ 2 1 z* ^ 


(e - e ) exp ( 




(A.5) 


We next substitute (53) and (57) into (A, 2) to (A.5) , In 

the integrands in (A. 2) to (A.5) f since J^(Zj^, z^) , ^2^ ' 

(Zj^t ^ 2 ^ and J^ (z^, Z 2 ) are slowly varying on the £ scale, 
1 2 

v/e replace them by I^(z^) , I^(z^) , I^(z^) , and I^*(z^) re- 

spectively. 

Because of tlie range restriction on z and z', most of the 
contributions of the integrals in (A. 2) to (A.5) come from terms 
with constructive interference. We conclude then that I,(z, z') 
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and small compared with ^^^ 2 , z') and I^(z, z’) . 

Balancing terms on both sides of equation (A.l', by their 
roepectivo phases, we arrive at the following equations 

u z 

ih,' (z - 2 ') irio^ - irio*3' 

Jy( 2 , z') e = e ^ l^izr z') 

“(ni + + ini2 - in,*z' 

^12 ® I^(z, z') A. 6 a) 




z*) e 


- 2*) “iriTZ + iTii*z 


1“ ’ -"‘‘I “1 
01 -"1 


~r«^e l_(z,z*} 


•f e 




I^(s, z') 


{A. 6b) 


Q 


+ =•) 


cr {2f s'5 ^ 


2») 




in^z -r 


I-, {z, s') 


2M 


•f R 


12 


e 


14(2, z») 


(A .6c] 



(z *^ z') 


J„ {z, z') e 




4- z M 


-ItItZ “ iri9*z' 

Rqi e I^^Cz, ZM 


-i(Ti * 4* Ti^*)d - irii^z* - iri„z 

H 2 ® ^ l4(z,z'). 


(A »6d! 


Setting z ' -»- z in equation (A,6) , we obtain equations (62) to 
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Because the fluctuating part of the permittivity is real, we 
can obtain from equation (6) the energy conservation relationship 

V 5 ^ ”|EJ^ =0. (A.7) 

We use superscript double prine to denote imaginary part of a 
quantity . 

Taking the ensemble average of (A.7) , we have 




For our purpose of energy conservation to the zeroth order, 
we require that both and y 2 to be of order 0(5) , We re- 

mark that y^ being of 0(S) does not imply y 2 is of 0(6) . 

For instance, the integration of y^ from z - -«» to 0, which 
is the case for a half -space medium, may yield a significant value 
for y 2 . We thus require that there be no constructive interfer- 
ence terms of 0(6) in y^. This condition then is sufficient 
to guarantee that y^ is of 0(6) , 

If 

Prom (52) , we obtain 


dz 




E* — E - E 

^^2 Hm ^Im 


dz' 


E* 

Im 


= -2i(ky 


E 


Im* 


+ 2i 


I 


— d 




(A.IO) 


In view of (55) and (17) , 


lim 

z-t-z* dz^ 


<Siz) <f*(z')> - <^*iz) (^(Z')> 




^^1 ^llm^^^ ^1^ Bff2, 2^) 
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4 . 2i 


.0 

—a 


(A. 11) 


Substitutdng Ca, 10) into (A.ll) , we find that = 0(6^)* Thus 
energy is conserved to the required order. 

If we had used the bilocal approximation to Dyson’s equation, 
rather than the nonlinear approximation, i ,e , instead of (13), we 
use 


Q(Z, Z^) ^ B^(2, 2^) 


<A,12) 


where ^1^ Green's function for the unperturbed 

problem, we would have obtained 




>= 2i 






(A. 13) 


We pse superscript b to denote bilocal quantities . If for 
2 - £ 0 (il) , 




(A. 14) 
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then as given by (A, 13) is of 0(6^), Equation (A. 14) is in 

general true for unbounded medium but not for bounded medium. We 
also remark that if we use bilocal quantities in the derivation, 
the MRT equations so obtained, even for half space medium, lack 
symmetry and do not conserve energy. 
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Figure Captions 


Figure 1 . 


Figure 2 , 


Figure 3 . 


Geoinetrical configuration of the problem, 

Einissivity of a half-space medium with = 3,2eQ^ 

= ,16eQ, % — 2 mm, 

Emissivity of a two layer random medium with = 3,2 £q, 
ej|^ = .ICSq, £ = 2 mm., ^ ~ 20cm. 
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APPENDIX VI 


Microwave Thermal Emission* from a Stratified Medium 
with Nonuni form Temperature Distribution * 

by 

L, Tsang, S. MjokU/ and J. A. Kong 

Department of Electrical Engineering and 
Computer Science and Research Laboratory 
of Electronics 
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Cambridge, Massachusetts 02139 


Abstract : 

Using the dissipation-fluctuation approach, the brightness 
temperature of a stratified medium with inhomogeneous permittiv- 
ities and nonuniform temperature profiles are solved. The solu- 
tion is expressed in closed form ready for machine computation. 
Numerical results are illustrated and compared with closed form 
analytical solutions and results obtained from the WKB method 
for simple profiles. 


* This work was supported by the NASA Contract 953524 and the 
Joint Service Electronics Program Contract DAAB07-74-C-0630 . 


I , Introduction 


In microwave passive remote sensing of the earth, the bright- 
ness temperature is determined by the geometrical configuration, 
the medium properties, and the temperature distributions of the 
area under observation. With the model of a vertically structure 

ril [21 

meditim Stogryn ^ used the fluctuation-dissipation theorem’- and 

f02miulated the solution for the brightness temperature in terms 

of a two-point boundary value problem by solving a second-order ■ - • 

ordinairy differential equation, together with the evaluation of 

an integral. In this paper we extend the propagation matrix for- 
[31 

mulation and solve for the brightness temperature of a stratified 
medium, The problem is an -important .one because inhomogeneous 
permittivities and nonuniform temperature profiles can be conven- 
iently approximated by a stratified medium. On the other hand 
analytical results can only be obtained for very few media with 
simple profiles. The solutions for the stratified medium are ex- 
pressed in closed forms and numerical examples are given for dif- . 
ferent temperature and permittivity profiles. The results are 
compared with the well-known exact solutions for simple profiles 
and with those obtained with the WKB method. 
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I3:« Dyadic Green *s Functions for a Stratified Medium 


Consider a stratified medium consisting of N layers and 
with planes of stratification perpendicular to the z-axis (Figure 
1) . Let there be a point source in region 0, In order to calculate 
field vectors in any region Jt, we first consider dyadic Green's 
function G(r, r') in the absence of the stratified medium. 


ECr, w) =Jd^r* S(r, r') • J(r'), 


( 1 ) 


The dyadic Green's function G(r, r') is governed by the equation 


. V 3c V X G(r, r') - G(r, r') = iwp I 6 (r - r') . (2) 

By the method of scattering superposition ^ , we obtain (from Ap- 
pendix A where dyadic Green's functions are derived for uniaxial media) 


Cf, r*) = 4 


^ CO 

“ ff 


dk..dk.. ^ {e(k^)e{k^) + h(k^)h(k^)}e 


iJ^ • (r - r') 


X y 


z ' < z 


ff dk dk -L {e(-j: )e(-k^) + h (-k,)h (-k,) }e^'' ’ 
Btt'' •»« ^ ^ k z z z z 


top 

p2 




z < z 
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(3) 
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where 


= (co^ye - 

e{k ) = -— - vyo " y 

^ (V + ky2)V2 y X 


h(k^) = i {e X k} ' ^ '■" ■ ' ' : ^ 

k . . ... 

* 

k = {k^^ + ky2 + k^2jl/2 

_ A n. 

• k == X k^ + y ky + z k^ . £ 8 ) 

„ A A A 

K = X k^- + y ky - z k^. (9) 

^ A ..... ^ 

We see that e Ck^) is a unit vector perpendicular to the plane 

— ^ ^ 

of incidence determined by k and z and h is a unit vector 
parallel to the plane of incidence . 

For the stratified itiediumr dyadic Green 's functions are 


written as 




_ _SE rr” 

8 ir^ 


dk dk — 
X y j. 


Ur'® eCVei’^ 


+ e(-k_)e^'' ■ *^)e(-k )e"^" ' + (R™ h(k Je^*' ' ^ 

Z 2 » 


A 

-f- h{-k )e 
z 


Ik 


^)h(-k^)e”^‘^ 



ao) 


for 2 < z ' in region 0 , 




X') » 


- ff 

8 ir^ 


ik, 


jj 

— « 


dk^dky 


((A, 






1 . 


A itCo • r ^ , ", 

+ .Djjhj, (-kj^^)e )h(-kjj)e } 


( 11 ) 


in region Z , and 




in region .t,. In writing the dyadic Green's functions, we use the 
first subscript to denote the region in which the pbservation point 
is located and the second subscript to denote the region in which 


the source is located. 


The coefficients C^, D^, and 


the reflection and transmission coefficients for ’TE and TM waves 
TE TM TE TM 

R , R , T and T are determined by matching boundary con- 

f31 

ditions and related to one another by propagation matrices , In 


this case, the TE and TM quantities are both for electric field 
vectors and are not dual to each other, ’ The TE quantities are 
found in reference 3, The TM quantities are given in Appendix B, 


-W-,' 


• * ■ Mv' 
» ..7jj I?** 
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III. Brightness Temperature of the Stratified Medium 


To find the brightness temperature of the stratified medium, 

E2] 


v?e apply the fluctuation -dissipation theorem . The thermal 
sources in region % is denoted by to) . The expection value 


of <Jj^ 


<Jjj^(r, oj) J^*(r', m')> = ~ coe^^" KT^(z)i6 (tn - w')5(r - r') (13) 

7T 


, f 


where K is the Boltzman constant, the temperature distri- 

bution in the £th region, and T the identity matrix. Equation 
(13) is true under the assumption hw << KT which is valid at 
microwave frequencies , 

The brightness temperature at angular frequency w is defined 


as 


[13 


K 


T^^Cic, w) = (c/u)^ “ f 


^00 


i(ic - ic') • r - i (oj - oj ') ti 


{e • <E(k, a))E*(k', aj')> • e e' 
for the horizontal polarization or TE waves and 


(14a) 


Li) ~ — (c/u)^ — ce [ dbi' f k^dk y d^k 

® K 2 0 Q 


{h • <E{k, to)E*(k', w')> • h e^^^ ” ^ *■ 


(14b) 


for the vertical polarization or TM waves. In order to determine 
the electric field intensity E(k, w) , we note that 
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.Od 


' t 

E(r/ w) - £ A dx 

Z = X -« 


r dy/"* - 1 


■•<09 


dz* G^^(r, r') . Jjj(?') (15) 


where t ~ n + 1 and d^ By the symmetric properties of dya- 

die Greenes functions, G^j^(r, r‘) = G^^ (r', r) , where the super- 
script t denotes transpose. In view of (11), we find 


oi 


irr r*) = - r d^k 5(k^ - - kj - k M 

Sr® -'-CO 2 X y 


ik • r 


2^ e 

y k. 


^ A ”iK- • r* A -ikj, • r' 

{e (k^) [A^e^ (-kj^^) e + (k^^) e 3 


“iK« • f* -ik„ • r' 

+ h(kjj) [Cjh^(-fcjj^)e * 


A Cfhaiige of variable of integration from Jc , k to -k and 

” . -T X 

-ky is made in arriving at (16) . Introducing (16) in (15) , we 
obtain 


E(k, w) = - ^ ^ ^ II «i3i'dy’ I ^ ^ dz' 5(k^ - /w^pe - k^^ - kj) 


8ir^ £ = 1 


X y 


, yr y. -iK- • T* yv -ik_ • T' 

— {e(k^) (A^ej^(-k^^)e t B^e^(kj^^)e 1 + h(k^) 


A -XK- • f' A -ikp • r' ^ _ 

t'=S,h;,(-kj,)e +D^hj(k^^)e ]}.J^(r-). (17) 
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Using (13) and (14) we find the brightness temperature to be 


til) = fc^oos 0 £ f ^ ^ dz' T, (z') 

^ = 1 


1 ^ “ik-^z* 

• • ■*» 
z 


( 18 a) 


T™(lc, til) = Ic’cos e £ [ 

ft ^ ^ 


t E« " r-d. _ ^ 

^ ^ dz* T (z>) 

^ = 1 


1 

\ 


In the derivation of (18) , we use the relation 6 [k - /o^ue - k ^ - k 

z X 

e= U(k ) COS 0 6 (k - oj/jlF) , where U(k ) = 1 for. k > 0 and 

A . Z 2 

U(k„) = 0 for k < 0, 
z z 

The calculation of the integral in (18) is straightforv/ard. 
Without loss of generality, we let the plane of observation be the 
X ”• z -plane and set k = 0. The temperature T. (z') is a con- 

y ^ 

stant in region Jl, (the case of non-constant temperature can be 
generalized easily) we find that for horizontal polarization 
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similarly! for vertical polarization. 
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''£2 1 


(C^e 


k 

cos 6 


- k. 


X 


+ k. 


ik« _d. 


) ( V > _ '“4 


■ 


- 1 >, 






'*' ’^x ^’’^t IJTM12 e 


^ol’^tl" 


2*=tz" 


( 20 : 


For homogeneous half-space with constant temperature, we 





T™ = U - r™)T^. 


( 22 ) 


In the derivation of (21) and (22) , we made use of the identities 

V = ’"tz + ^t" ’'tz'/Hz- 

In Figure’ 2 we examine the angular dependence of the brightness 

TE TM 

temperature. At 1 GHz, Tg and Tg are plotted as a function 
of observation angle for the following profiles: 


e, tz) 


- 9,0(1 + i 0.3) - (5.5 + i 2.5) e 


az- 


(23a) 




(23b) 


T(z) = + AT e 


bz 


(23c) 


= T(-d) 


(23d) 


with a = 0,02 cm”^, b == 0.05 cm“^, d = 30 cm, and AT = ± 20®K, 

It is interesting to observe .a. maximum for vertical polarizations 
similar to the Brev/ster angle for a unifom half space medium. 

In Figure 3, we plot the brightness temperature of a radio- 
meter looking frcm nadir as a function of frequency for the following 
permittivity and temperature profiles: 


az 

-i = 9,0(1 + i 0.3) - (5.5 + i 2.5) 


(24a) 




<24b) 
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(24c) 


T - T + AT 
o 

T^ - T(-d) - (24d) 

with a = 0.02 cm“^, d = 30 cm', = 280°K, AT + 20®K, and 
b s= 0«05 cm"^ and 0.1 cm"^. In the calculations the temperature 
and the permittivity profiles are stratified into 350 layers 
from 0 to 2 = -d. We see that at low frequencies the subsurface 
temperature affects T^ more than at high frequencies , At very 
low frequencies the brightness temperatures for different parameters 
of a and b approach the sam© value dictated by T^ and the 
effective emissiyity of the medium. 

In Figure 4 the brightness temperature as observed from nadir 
is plotted as a function of frequency for the following profiles: 


e^iz)/Zo ^ 9.0(1 + i 0.1) ~ (5.5 + i 0.83) (25a) 

e^(z -20 cm)/e^ ~ 9.0(1 + i 0 .05) (25b) 

T{ 2 ) = T + AT e^^ (25c) 

o 

T(z = -20 cm) (2bd) 


V7ith a - 0,02 cm“^, b ~ 0.1 cm"^, = 280°K, and AT = ± 20®K. 

Note that there is an abrupt change in the permittivity profile 
at z = -20 cm. At low frequencies, we see interference effects. 


For higher frequencies, the brightness temperature is determined 
by ± AT times the effective surface emissivity of the medium. 

• The interference effect at low frequencies is not character- 
istic of the abrupt change of the permittivity profile. For a 
continuous permittivity profile the interference is illustrated 
in Figure 5, where 

ei(z)/EQ = (2.88 + i 0.34) (26a) 

= Ej^(z - -1 m) (26b) 

T(z) = 300 - 20 -CO < z < 0, (26c) 

“1 -1 

a = 2 m and b = 3 m , In the calculations, the region from 

z = 0 to z - -1 m is stratified into 350 layers. The tempera 
ture profile from z = -1 m to -® is not stratified because 
(18) can be readily integrated. For, this profile the exact solu- 
tions are presented in Appendix C. The result calculated with the 
stratified model is found to agree very well . with that .calculated 
with the exact formulas . 
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IV. WKB Solution 


When the permittivity of the stratified medium can be written 
in terms of a continuous function e(z) for z < 0, we foJ.low 
an approach parallel to that used by Stogryn and find the brightness 
temperature for the horizontal and vertical polarizations to be 

T™ = — / dz T(z) e “(z)|i|i(z)i=' C27) 

and 


tTM f dz T(z) z"iz) 

.B V ^ 


where 


1 d[J) {%) 

2 


2 

e^(z) dz 

IT 




(2 8 ) 


G^(z) = (z) 




(29) 


is the dielectric constant. The functions ij;(z) and (f>(z) in 
(27) and (28) are governed by the following equations: 


— + [e^(z)k^ -k 23^(35 ^0 

dz® ^ ^ 


(30 a) 


E^(z) 


dz 


1 d(^(z) 


e^(z) dz 


+ [Ej,(z)k® - k^®]({)(z) = 0 


. (31a) 


The boundary conditions are 


1 + = ^( 0 ) 


(30b) 
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z - 0 


(30 c) 


k^(l - - 1 


3t}; 


1 + = 4-CO) 


E (0) k (1 - RJ 
r . z V' 




Approximate solutions to these equations are easily 
with the WKB method. 

Consider an equation of the following form 


+ [g(z) l^X = 0. 

dz^ 


!Ehe WKB solution for an outgoing wave is 


J dz g(z) 

Xtz) — ■ e 


Provided 


g » i ^ » Jl (l//g) , 

- g dz »/g dz^ 


Comparing (30) and (31) with (32) , we see that 


+ Ig„ (z) ]^y = 0 

dz=^ ® 


tg^(z) 1^ “ 


reproducibility of the 
ORICjINAL page is poor 


and 


(31b) 

(31c) 

determined 

(32) 

I 

(33) 

(34) 

(35) 

(36) 
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+ [g_(z) ] 2 [«f)//e“] = 0 ( 37 ) 

dz^ ^ ^ r 

Egj„(z} 3 ^ - Is^iz)]^ - (dVdzMa/ZH^). ( 38 ) 

The condition ( 24 ) is met if the medium is slowly varying such 
that 

— — C/^) « g^(z) . (39) 

dz ^ ® 

r 

Under this assumption, 

[gjjj(z)]® ^ [g^Cz)]® = ( 40 ) 


After matching boimdary conditions at z = 0 ar;d making use of 
( 27 ) and (2 8 ) , we obtain 


= 

B 


(1 - r 


TE 


) / dz 2 T(s) SlilL g«( 2 ) 

I g(z) i g' (0) 


2 f dsg"(s) 
0 


( 41 ) 


mmCO 


E 


(1 - r™) f dz 2T(z) 


2 -f dsg"{s) 

f (z) |g( 0 ) [ g' (z) „ > . e ' ° 

f( 0 > |g(z,| g.( 0 ) " ‘ ’ ■ 


( 42 ) 


where 
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g 




kg - g(0) 
•}■ gCO) 



k^ Ej.(0) - g(0) 
kg Ej,(0) -l- g{0) 


2 


(4 


jjj. 


(43c) 


It is interesting to see that under the assumption that g" is 
very small such that |g(z)| g'(z)/ (41) and (42) reduce to 


• T. 


B 


(1 - r) / 2T(z) g"( 2 ) e ° 


r ^ 

2 / dsg”(s) 

J r\ 


(44) 


which is a well-known result in the radiative transfer theory. 

To illustrate the use of (41) and (42) , we plot the V7KB 
result in Figure 5 to compare with that obtained from the stra- 
tified model. We note that the WKB result is quite accurate 
in the high frequency side. For low frequencies, the V?KB fails 
to predict interference effects. This is due to the fact that 
the reflectivity is a constant for the WKB approach, whereas 
for the stratified model, the reflectivity exhibits the inter- 
ference effect. The reflectivities are plotted in Figure 6, The 
CPU time in using the IBM 370 for calculation of the brightness 
ten^erature in Figure 5 is approximately 0.12 minutes for the 
stratified model, 0.10 minutes for the exact solution, and 0.03 
minutes for the T\TCB approach. This implies that t^nKB is useful 
in the frequency range where its results are valid. 






V. Discussions 


In, this paper we apply the dissipation-fluctuation approach 
to find the brightness temperature of a stratified medium with 
inhomogeneous permittivities and iionuniform temperature profiles. 
The results are expressed in closed form. The solution to this 
problem is an important one because only very few special profiles 
can be approximated and solved in analytical forms. Even then 
the solution is in the form of special functions and requires per- 
haps a.9 much computer time to calculate numerically. Numerical 
results to the stratified medium model are illustrated. They are 
also compared with the exact solutions and the results obtained 
from WKB approximation. At high frequency limit, profiles that 
can be treated X'/ith the WKB method may require less computer 

time. However when the frequency gets higher than microwave f re- 

15 1 

quency, the approach proposed in this paper becomes inaccurate 
Also in the calculations we omit the contributions to the bright- 
ness temperature arised from reflections of sky radiation and 
cosmic background, which can be easily accounted for by adding to 
the above results . 
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Figure Captions 


Figure 1 . Geometrical configuration of stratified isotropic 
medium. 

Figure 2 . Brightness temperature as a function of observation 
angle for and TjM waves at 1 GHz, 

Figure 3 . Brightness temperature as a function of frequency 
for profiles in Equation (2 4) . 

Figure 4 , Brightness temperature as a function of frequency 
for profiles in Equation (25) . 

Figure 5 , Brightness temperature calculated with the stratified 
model and the exact formula compared with that cal- 
culated v/ith the WKB method. 

Figure 6 , Reflectivities calcuated for the stratified medium 
and for the WKB method. 
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Appendix A 


In this section we derive dyadic Green’s function for a 
uniaxial medium with optic axis parallel to the z axis , Dyadic 
Green’s functions for uniaxial media are governed by the following 
equation 



G(r, r') = iw 16 {r ~ r') 


(A.l) 


where 



1/V. 0 0 

0 x/u O' 

0 0 l/y^ 


is the inverse permeability tensor. 


(A. 2) 



{A. 3) 


is the permittivity tensor, I denotes the identity matrix, and 
V is defined by V • A = V x a for any vector A. Equation (A.l) 
follows directly from Maxv;ell's equations. 

We write the Fourier transform pair of dyadic Green's function 
as 
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(A. 4) 


S(k, r*) -/ d^r 6(r, r’) * 


09 


5(r, r') = f d^k G(k, r') 

&TT ® "® 


(A .5) 


Substituting (A, 5) in (A.l) yields 


G (k , f * ) = -iw [K 



k + w^e]~^ 



r • 


(A. 6} 


vjhere 


k 




(A. 7) 


After introducing (A.2) , (A.3) , and (A. 7) into (A.6) and performing 
the matrix inversion, we obtain 
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where 


D^(k) 



e) 


CA.9) 


~ (k^® -f- ky^ - to^ye^) . 


{A. 10) 


To evaluate the integration over k in (A. 5) , we find that the 
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poles of the integrand occur at the zeros of D 2 <k)f 

corresponding to = ± k^^ and k^ *= ± k^® where 


m 


kj,- = - y(kjj^ + ky=>)/p. 


(A. 11) 


and 


k^® = - e (k^^ + k ^ • 


(A. 12) 


[31 

They correspond to the type I and type II characteristic waves 
in the Jnedium, Assiming the mediiim is slightly lossy, then we find 


- J2E /■/" 


(r, r ') = •{ 




1 

ee e 


m 


+ AhO'z )h(k/)e 

_ 

for z ' < z 


f^//“ 


8tt 




^ ee e 


(r “ r') 


m 


L 


1 1= e> 

+ _1_ h(-kj,®)h(-kjj®)e 


(r - r') 


for z < z ‘ 


where 


(A. 13) 
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(A. 14) 


= X 4 y ky + z 
— ^ ^ in 

* ’'x y ’'y * ’'z 

/V A A 0 

Kg = k kj. 4 y fcy - z kg 

K - K k + y k - z k ^ 
m X -* y z 


{A, 15) 


(A. 16) 


(A.17) 


A 

e 



(X ky - y kj^) 


h(k ®) =z 1 [e X k 3 

^ k ® 


{A. 18) 


(A. 19) 


h(“k ®) rz 1 Ee X K ] 
z j, e 


CA.20) 


and k = u/pe*. We see that e is a unit vector in the direction 

A 

of the electric field for TE waves and h is in the direction 
of the electric field for TM waves. Equation (A. 13) reduces to 
Eq. (3) when the mediuin is isotropic. 





The TM reflection coefficient in continued fraction form is 


where 


, \ -i2k,.d, . 


, \ -i2k,_ 


(d, - d,) 


^ Ch - i)n 


^nz' n n ~ 1 


<n - l)n 


+ 1 ) 


= + 1 ~ t^l) z 

+ 1 ^(A + l)z 


is the TM Fresnel reflection coefficient at each boundary. 
The wave amplitudes and in each layer are related by 

propagation matrices . 
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Appendix C . 

In this; appendix, we present the exact solutions of the 
brightness temperature of the medium with permittivity profile 

£3^(2) = (£3^* + i 

e^Cz) = ~ -d) 


and temperature profile 
T(z) = + (AT)e^^. 


In terms of Bessel functions, the brightness temperature as 
observed from nadir is equal to 



+ T 


It * 
''t 


+ s Hotna) ,P 


.2k." 

+ t f^\ 

(b + 2k") It ] 

- 1 ' o 


-bd 


where 


Trrio ^ 
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